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ABSTRACT 

Aims. We propose a new method of pushing Herschel to its faintest detection limits using universal trends in the redshift evolution of 
the far infrared over 24 pm colours in the well-sampled GOODS-North field. An extension to other fields with less multi-wavelength 
information is presented. This method is applied here to raise the contribution of individually detected Herschel sources to the cosmic 
infrared background (CIRB) by a factor 5 close to its peak at 250 pm and more than 3 in the 350 and 500 pm bands. 

Methods. We produce realistic mock Herschel images of the deep PACS and SPIRE images of the GOODS -North field from the 
GOODS -Herschel Key Program and use them to quantify the confusion noise at the position of individual sources, i.e. , estimate a 
“local confusion noise”. Two methods are used to identify sources with reliable photometric accuracy extracted using 24 pm prior 
positions. The clean index (Cl), previously defined but validated here with simulations, which measures the presence of bright 24 pm 
neighbours and the photometric accuracy index (PAI) directly extracted from the mock Herschel images. 

Results. Both methods converge to comparable depths and fractions of the CIRB resolved into sources individually detected with 
Herschel. After correction for completeness, thanks to our mock Herschel images, individually detected sources make up as much as 
54% and 60% of the CIRB in the PACS bands down to 1.1 mjy at 100 pm and 2.2 mjy at 160 pm and 55, 33, and 13% of the CIRB 
in the SPIRE bands down to 2.5, 5, and 9 mjy at 250 pm, 350 pm, and 500 pm, respectively. The latter depths improve the detection 
limits of Herschel by factors of 5 at 250 pm, and 3 at 350 pm and 500 pm as compared to the standard confusion limit. Interestingly, 
the dominant contributors to the CIRB in all Herschel bands appear to be distant siblings of the Milky Way (z~0.96 for T<300/rm) 
with a stellar mass of M*~9x 10 10 M 0 . 

Key words, infrared: diffuse background - infrared: galaxies - galaxies: statistics - galaxies: photometry 


1. Introduction 


Primeval galaxies in the process of forming their first gener¬ 
ation of stars were predicted to generate a diffuse extragalac- 
tic bac kground light (EBL) in the optical and near-infrared 
(NIR) ( Partridge & Peebleslll967l). la ter on extended to the far- 
IR (iLow & Tucker 11*9681 Kaufmanlll976h . iPuget et aD (l996) 
reported the first measurement of the cosmic IR background 
(CIRB), an isotropic IR signal appearing in the all-sky maps 
of the FIRAS experiment onboard COBE with an integrated in¬ 
tensity of ~30±5 nWm 2 sr 1 between 100 and 1000 pm. The 
CIRB represe nts about the half of the EBL dDole et al.1 l2006t 
lLagache et alJ 120051) . Although it sums up only ~2.3% of the 


intensity of the cosmic microwave background (CMB), it is the 
strongest emission due to galaxies so that it represents the net en¬ 
ergetic budget of galaxies. The cosmic UV-optical (COB) back¬ 
ground and CIRB are expected to be generated by the energy 
produced through stellar nucleosynthesis and, for a minor frac¬ 
tion, by black-hole accretion and emission from the diffuse in- 
tergalactic medium (lHauser & Dwekll2001l) . The contribution of 
accretion in the form of active galactic nuclei emission is ex¬ 
pected to represent a fraction of ~ 15% of the CIRB dJauzac et al.l 
1201 lb . Since the discovery of the CIRB, a major issue in obser¬ 
vational cosmology has been to identify the individual sources 
responsible for this diffuse light which contains the net energetic 
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budget of nucleosynthesis integrated over cosmic time, hence the 
detailed information on how and when stars formed in galaxies. 

Several factors make the correct identification of sources a 
challenge when multi-wavelength data of extragalactic objects 
are combined. The resolution 0 of a telescope having a collect¬ 
ing mirror of a diameter D is limited by diffraction, i.e. 6 ~ A/D, 
a problem particularly difficult to treat for far-infrared (FIR) and 
submillimetre instrumentation (A > 100 gm). It is possible to 
attack this problem by constructing larger collecting surfaces 
for the telescopes (as well as instruments with smaller detec¬ 
tion elements) until the desired angular resolution is reached. 
However, this increases substantially the associated costs that 
scale with the size of the telescope mirror to the 2-2.5 power 
(Ivan Belle et al.ll2004b . as well as rocket limitations that prevent 
the launch of big telescopes in space. Poor spatial resolution 
plays a major role in confusion noise, which can be defined as 
the uncertainty in the measurement of the flux of a point source 
owing to the variance of the underlying background induced 
by extremely faint sources or by the contamination of nearby 
point sources. The origin of the background may be instrumental 
(caused by the thermal radiation of the telescope itself), intrin¬ 
sically dif fuse (caused by the inter stellar medium, e.g. galactic 
cirrus, see iHelou & Beichmanil 19901) or may be the result of the 
accumulated light of faint objects (the extragalactic background 
light). 

The launch of the Herschel satellite dPilbratt et al.ll2010l) has 
opened a new window on the FIR Universe by allowing the 
detection of individual galaxies down to uniquely faint limits 
with deep extragalactic surveys (for a comprehensive review see 
iLutzl 120 14h . Thanks to its 100-500 gm coverage, Herschel pro¬ 
vides the best insight into the sources responsible for the peak 
emission of the CIRB. Although it is the largest telescope ever 
sent into space, its 3.5m-class observatory does suffer from the 
so-called confusion limit, and the depth of these extragalactic 
surveys is essentially limited by the fluctuations of the signal 
due to a sea of unresolve d sources, together with the blend¬ 
ing from bright neig hbours dBerta et alJ2010t lOliver et al.l2010t 
Ngu ven et al.l l2010). Confusion is inherent to IR data, particu¬ 
larly in the 250-500//m SPI RE bands, where the confusion limi t 
is as high as ~20 mJy dOliver et alJl2010t iNguven et al.ll2010h . 
but in this work we attempt to develop a method to go beyond 
the nominal confusion limit. The confusion noise results from a 
combination of effects including the beam size, the instrumental 
noise, astrometric uncertainties, and the emission of faint, unde¬ 
tected sources. Isolating the effect of each of these contributions 
to the confusion is a challenging task (e.g. degeneracy between 
the position of the sources and the effect of the uncertainty in 
individual pixels produced by the instrumental noise, pixel sen¬ 
sitivity, and local fluctuations due to faint sources). 

From this point, this paper is divided into the following sec¬ 
tions. Section [2] explains the new method to incrementally cor¬ 
rect the observed Herschel colours from shorter to longer wave¬ 
lengths; §[3] presents the ancillary dataset supporting this work; 
$4] analyses the role of clean galaxies, blending and clustering 
effects in the evolution of mid-IR (MIR) to FIR colours as a 
function of redshift; §0 presents the relations that best fit the 
colour-redshift relations, as well as the iterative corrections ap¬ 
plied to correct blending effects; §j6] describes how the realistic 
mock images were made; §0 provides an analysis of the con¬ 
fusion noise in the extracted sources from the simulated maps. 
Finally, ^8] presents new values for the resolved CIRB and the 
typical redshifts and stellar masses of the galaxies responsible 
for its origin, and §|9]summarises the main results of this work. 


2. The method 


This work relies on the combination of the information of Spitzer 
and Herschel from the shortest to the largest wavelengths thanks 
to the existence of universal mid-to-far infrared colours of galax¬ 
ies as a function of redshift. The success of this propagation of 
the information from the shortest to the longest wavelengths re¬ 
lies on the recently found evidence of universality of the IR spec¬ 
tral energy distribution (SED) of galaxies (iElbaz et al.ll201lf) . it¬ 
self reflecting more general scaling laws now often referred to 
as a “main sequence of star-forming galaxies”. Recent evidence 
suggests that most galaxies grow in stellar mass at a steady- 
state rather than from episodic, merger-driven starbursts, and 
new scaling laws have emerged from this scenario of galaxy for¬ 
mation. 

First, a tight correlation between the specific SFR ( sS FR - 
S FR/M+) and the stellar mass M+ at different re dshifts 
z ~0 dBrinchmann et ah 20041; ISalim et al.1 12007). z~l 


(Elbazetal||20Q7^ Noeske et al .1120071). z~2 (Daddi et al. 2007b; 


Pann ella et al.l l2009l). z~3 ' (Magdis et al.1 120 1 0t Elbaz et al.1 


201T), z~4 ( Daddi et al.ll2009l) . In this scaling law, most galaxies 


fall on a main sequence (MS) with steady star formation, while 
a small fraction of outliers experience starburst (SB) episodes. 

Second, there is a correlation between the total IR luminos¬ 
ity L/s and the luminosity at 8//m, L X; , m . The L/« is mostly due 
to the FIR emission of big grains of dust, while the 8pm emis¬ 
sion is dominated by polycyclic aromatic hydrocarbon (PAH) 
molecules. Again, most galaxies fall on a main sequence with 
a median //(8 = T/ W /L X;lm of ~4.9. Outliers to this sequence, 
i.e. T/«/T8/jm >>4.9, are associated with mergers that show com¬ 
pact star-forming regions, where PAH molecules are more effi¬ 
ciently destroyed hence decreasing in L 8/Jrn dElbaz et al.ll201 It) . 

The third scaling law separating a dominant mode of star for¬ 
mation from a minor mode of starburs ts is provided by the ex¬ 
tension of the Schmidt-Kennicutt law (ISchmidi 19591; IKennicutt] 
1 19981) . the tight relati on between projected densities of SFR and 
gas, to high redshifts (iDaddi et al'll2010HGenzel et al.ll201 (H) . 

As a result, the mode in which galaxies form their stars 
seems to follow a set of simple rules. The infrared SED of dis¬ 
tant galaxies is on average very similar to the one of local galax¬ 
ies and these simple scaling laws that have been found can be 
used to extrapolate total IR luminosities from single photomet¬ 
ric measurements. In this work, considering the potential of ex¬ 
trapolating the IR fluxes from shorter wavelengths and taking 
advantage of the unique multi-wavelength data existing in the 
GOODS -North field, we present a method of digging below the 
standard global confusion limit (i.e. the statistical limit that does 
not account for the actual positions of individual sources) and 
reach a local confusion limit, accounting for these positions, to 
resolve most of the cosmic IR background even in the SPIRE 
bands that are most affected by confusion. 

We achieve this by building realistic simulated images of a 
cosmological deep field that is well-sampled in redshift and has 
a very complete multi-wavelength dataset. The robustness of the 
source extraction and the accuracy of the photometry are studied 
in the simulations, aiming to find a criterion that determines the 
reliability of fluxes by 


- quantifying the impact of confusion on the uncertainty in the 
determination of Herschel flux densities 

- understanding the relative impact of the different sources of 
noise 

- improving the selection of prior sets of target positions to im¬ 
prove the search for faint galaxies in the process of catalogue 
building. 
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The resulting criteria are designed to be applicable to other cos¬ 
mological fields that may lack a complete and well-controlled 
redshift, and multi-wavelength ancillary datasets. 

Consequently, we proceed according to the following steps 

- We exploit the fact that prior sources, seen at shorter wave¬ 
lengths with higher spatial resolution (e.g. 24 pm sources), 
are located at positions that are tightly correlated with the 
presence of FIR sources. 

- Thus, instead of defining a global confusion noise (produced 
by a floor of faint sources percolating over the full image, and 
to limit detections above 3cr of this global noise) we use our 
simulations and the information at the shorter bands of each 
given Herschel band to attribute to each source a photomet¬ 
ric accuracy index. This index measures a “local confusion 
noise”, i.e. the impact of confusion at the actual position of 
individual sources, by estimating the difference between the 
input and measured flux densities of every individual source 
used to produce the mock Herschel image. 

- We use this index to isolate a subsample of sources that are 
the least affected by confusion from neighbouring sources. 

- We define this revised detection limit for a subset of the 
sources as the local confusion noise. 

- We are then able to dig into the heavily confused deep FIR 
images, improve the Herschel detection limits, resolve a 
greater fraction of the cosmic IR background (CIRB), and 
determine what fraction of the CIRB can be resolved into 
individual sources. Since sources that are rejected for being 
affected by random projections of sources in the line of sight 
are not biased, the properties of the sources with clean pho¬ 
tometry are expected to be representative of the whole sam¬ 
ple. 


3. The GOODS-North dataset 

The full 10' x 15' Great Observatories Origins Deep 


(Poglitsch et a 

1. 2010) 

(Griffin et al. 

2010) at 


at 100 and 160 pm, and SPIRE 


GOODS -Herschel Key Project (lElbaz et al.ll201 ll) . By construc¬ 
tion, the SPIRE images properly encompass the GOODS-N field 
but also cover a wider area that is included in the FITS images of 
the present data release but not in the catalogues due to the lack 
of ancillary data in the external parts, in particular deep 24/jm 
catalogues that are required to perform prior Herschel source ex¬ 
traction. The total observing time in GOODS-N is 124 hours and 
31.1 hours for PACS and SPIRE, respectively. 

The mosaics used in the present paper are the same ones as 
those described in lElbaz et al.l (1201 ll) . In short, the data were re¬ 
duced within the HIPE environment in the Herschel Common 
Science System (HCSS) with the official PACS Photometer 
pipeline dWienrecht et alj 120091) . Maps were created from the 
timelines of each AOR via the HCSS photProject algorithm, 
which is equivalent to the drizzle method (Fruchter & Hook 
2002). Each individual AOR map was projected with the same 
world coordinate system (WCS). Since the SPIRE bolometers do 
not fill the field of view, scanning at a particular angle is neces¬ 
sary in order to provide relatively homogeneous coverage, hence 
the observations were executed as cross-linked scan maps. Each 
single SPIRE cross-linked map covers 22' x 19' with 3 repeti¬ 
tions for a total field of view of 900 arcmin 1 2 3 . 

The pixel sizes are 1.2 arcsec and 2.4 arcsec for the PACS 
maps at 100 pm and 160 pm maps and 3.6, 4.8, and 7.2 arcsec 
for the SPIRE maps at 250, 350, and 500jum, respectively. These 


Spitzer 

A (pm) 

3.6 

4.5 24 

70 


FWHM 

1.6” 

1.6” 5.7” 

18” 


Herschel 

A (pm) 

100 

160 250 

350 

500 

FWHM 

6.7” 

11” 18.1” 

24.9” 

36.6” 


Table 1 . Full Width Half Maximum (FWHM) of the Point 
Spread Function (PSF) of the Spitzer IRAC and MIPS bands 
(upper part), and Herschel PACS and SPIRE bands (lower part). 


significantly oversample the PSF full width at half-maximum 
(FWHM) by a factors of 5, corresponding to a sampling finer 
than Nyquist by a factor 2 (see Table Q] for a summary of the 
FWHM of Spitzer and Herschel.) 

The present SPIRE data from the GOODS -Herschel project 
in GOODS-North are the deepest data ever obtained at 200- 
600 pm, together with those similarly deep data obtained by 
the HerMES '] consortium (Oliver et al]|2010l) in GOODS-South. 
We note that deeper data exist in the GOODS-South field for 
PACS for the central 8' x 8' region from GOODS-Her.se/je/and 
extended over the whole GOODS-South field in the PEF0- 
GOODS-Herschel combination dMagnelli et al.l[2013l) . The 3-cr 
depth for the PACS data in GOODS-North are 1.1 mJy (100 pm), 
2.7 mJy (160pm) as compared to 0.6 mJy (100pm), 2 mJy 
(160 pm) after combining the GOODS -Herschel data with the 
PEP ones dMagnelli et al.l [ 20 T 3 I PI Due to the confusion limit, 
the data at 160 pm are mildly deeper, while the 100 pm data 
in the central region of GOODS-South are nearly twice deeper. 
However this deeper level is obtained at the centre only of the 
field and at the expense of a varying coverage on the outskirts 
of the field. Since the present paper presents a method to deal 
with confusion, and that this effect only starts to play a dominant 
role with Herschel above 160 pm, the following results would 
be comparable if done in one or the other field. Forthcoming pa¬ 
pers will apply this method to all the CANDELS fields, i.e. the 
two GOODS fields, as well as the central part of COSMOS and 
the UDS field in the framework of the Astrodeep FP7 European 
project. 

3.1. Source extraction method 

Source extraction on the real Herschel images, as well as on 
the mock images that will be later on described was per¬ 
formed, as described in Section 4.1 of Magnelli et al. (2013, 
see also lMagnelli et akll2009l) . A MIPS-24pm catalogue (.S - 2 4;lm 
> 20 pjy) w as built from a list of so urces coming from a blind 
SExtractor dBertin & Arnoutsl 11996 ) source extraction on an 
IRAC mosaic made from the weighted average of the 3.6pm 
and 4.5 pm images, using a ‘Mexhat’ convolution kernel to en¬ 
hance deblending in crowded regions. Due to the pixel sizes of 
the IRAC camera, the spatial resolution of both IRAC bands is 
given by the same FWHM-1.6". 

As a result of this choice, the position uncertainty of the 
Herschel sources is the one of the IRAC priors but the accu¬ 
racy of the separation of neighbouring sources is controlled by 
the actual FWHM of the PSF at the Herschel wavelengths (listed 

1 Herschel Multi-tiered Extragalactic Survey 

2 PACS Evolutionary Probe (PEP) Guaranteed-time Survey 

3 As explained in Section 5 of lMagnelli et all ( 120131) . although the 3- 
cr detection limit is 1.3 mJy at 160 pm, sources below 2 mJy are affected 
by confusion, hence reliable individual detections are limited to a depth 
of 2 mJy rather than 1.3 mJy. 
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in Table Q}. If we had a perfect knowledge of the pointing of the 
satellite, and in the absence of instrumental noise, sources could 
be deblended with a nearly infinite accuracy, which we actually 
tested. But in the presence of these two sources of uncertainty, 
even accurate prior positions are limited by the fact that we are 
not certain of the position in the Herschel maps and due to the 
insttumental noise, which introduces an extra uncertainty on the 
centroid of the PSF that is fitted. The resulting flux uncertain¬ 
ties are therefore very difficult to quantify, requiring the most 
accurate mock Herschel mosaics to reproduce these sources of 
uncertainties to quantify precisely the difference between input 
and measured flux densities for sources distributed with the ac¬ 
tual source positions in the field. This is the key bottleneck of 
Herschel extragalactic surveys that we are addressing here. Note 
that stacking on a larger number of source positions does re¬ 
duce the effect of blending, since only the central sources are 
repeated, while neighbouring sources are diluted in the stacks, 
but that this is at the expense of studying large galaxy popula¬ 
tions with no insurance on the actual fraction of sources that are 
making the signal. We are here interested in digging deeper into 
the confusion regime but following individual sources that can 
later on be studied. 

We used MIPS-24gm priors from iMagnelli et all (12013l) im¬ 
posing a minimum flux density of 20 pJy (3-cr limit, while for 
PACS-160jt/m and SPIRE-250/im we restricted the 24 pm pri¬ 
ors to fluxes larger than 30/rJy (5-cr limit) (reducing the num¬ 
ber of priors by about 35%). For SPIRE-350 and 500 pm we 
used the 24 /mi prior positions for sources with a S/N < 2 de¬ 
tected at 250 pm. These criteria were adopted from Monte Carlo 
simulations to avoid using too many priors that would result in 
subdividing flux densities artificially, while producing residual 
maps (after PSF subtracting the sources brighter than the detec¬ 
tion limit) with no obvious sources remaining. 

It is important to avoid using a too large number of prior 
positions since this would lead to an “over-deblending issue’’, 
namely the flux density of a single Herschel source may be dis¬ 
tributed between several priors with high uncertainty. Therefore, 
in the event where two 24 pm sources would turn out to be closer 
than 1.5 pixel ( i.e. FWHM/3.5 or Nyquist/1.5), we only kept the 
brightest one for the source extraction. This choice does not im¬ 
pact the source extraction itself - it is only used to avoid over¬ 
estimating the noise in the measurement by trying to deblend 
too nearby sources - but it means that we have favoured one out 
of the two candidates as counterpart. The implications of this 
choice will be the matter of a forthcoming paper where we show 
that in a minor fraction of the cases, this approach leads to the 
missing of distant SPIRE sources that are wrongly attributed to 
more nearby galaxies (Shu et al. 2015, in prep.) 

The MIPS-24yum data probe fainter IR galaxies thanks to 
a less intense ^-correction resulting in a flatter redshift evolu¬ 
tion of the minimum L/« that can be detected. The Spitzer-MIPS 
24 pm image reach deeper confusion levels (20 pi y for 3<x) up 
to z ~ 3 than Herschel data in the GOODS fields and most 
Herschel sources have a 24pm counterpart dElbaz et al.|[201 ll) . 
A minor number of Herschel objects (21 galaxies, <1%) have 
no 24 pm counterparts. Those sources have been identified as 
MIPS-dropouts, and show strong silicate features in absorption 
at 9.7 and 18/um falling in the MIPS-24 pm band when the galax¬ 
ies are at z ~1.3 and z ~0.4, respectively (iMagdis et al.11201 lh . 
This means that source extraction based on 24 pm prior posi¬ 
tions can recover most of the sources in the Herschel images. 
It is expected that above a redshift of z ~3, the SPIRE bands 
will be more sensitive to distant sources at the depth reached 
here in GOODS -Herschel than MIPS-24 pm. As we will show in 


Section ROl this will concern only a small fraction of the sources 
that will not play a dominant role in the confusion limit. 

3.2. Flux uncertainties 

Two methods have been classically used in pre vious studies 
(e.g. iMagnelli et al.ll2009l 120131: lElbaz et 311201 lh to determine 
measurement uncertainties on the Herschel flux densities. The 
first one consists in measuring the noise level in a “residual map” 
at the position of the source. To produce such residual map, a 
first source extraction is performed to determine the 3-cr detec¬ 
tion threshold, then sources brighter than 3-cr are identified, and 
only those ones are kept in a second iteration of the source ex¬ 
traction code on the original images. This leads to the production 
of a residual map where only 3-cr sources have been subtracted, 
thus avoiding subtracting highly uncertain sources that would 
produce noise artificially. Keeping this residual map as a noise 
map, a third source extraction is then performed using all prior 
positions as in the first iteration but keeping the previous resid¬ 
ual map as a noise map. This noise measurement is our “local 
residual noise level”. 

A second method consists in injecting a small number of 
sources with a given flux density in the real Herschel mosaic, 
extracting the sources using their all the prior positions, and in¬ 
cluding the ones of the extra sample, for which positions have 
been blurred in the input prior list by a random uncertainty mim¬ 
icking the position uncertainty from the satellite. Reproducing 
this step many times for each flux density, we obtained an uncer¬ 
tainty on the measurement that we call “global simulated noise 
level”. This noise level can be compared to the local residual 
noise level for any of sources of similar flux densities. When 
this comparison was done, it was found that the two noise levels 
were equivalent for PACS, while for SPIRE the local residual 
noise level was underestimating the global simulated noise level 
by a factor 3/5, i.e. sources expected to be measured at the 5-cr 
level, were in fact detected at 3-cr. The inconvenient of this mea¬ 
surement is that it is not local, it is an average measurement over 
several positions on the map. 

In the present paper, we will introduce a better determination 
of the noise level that will be defined from the use of realistic 
Herschel mosaics. First, these mosaics will be used to demon¬ 
strate that the previous two noise levels can be made consistent. 
We will show - in Figures fl 21 13l discussed in Section 1731 - that 
the factor 3/5 is indeed reproduced in the realistic mock mosaics. 
Second, we will be able to determine a “local simulated noise 
level” by determining the noise that affects the measurement of 
a source at its real position but on the mock mosaic. This local 
simulated noise level will be determined by a “photometric accu¬ 
racy index” (PAI) that will compare the measured and input flux 
densities of sources in our simulated images, where all sources 
are injected rather than a small sample randomly injected at posi¬ 
tions that do not reproduce the actual local confusion noise level 
of sources. 

3.3. Redshifts and stellar masses 

Spectroscopic redshifts come from [Barger et al.l ( 2008 ) and 
Stern et al. (in preparation). 

Photometric redshifts were obtained using EAZ'tQ adopting 
a linear combination of the seven standard EAZY galaxy tem¬ 
plates as described in lPannella et al.l ( 2014 ). When compared to 
the spectroscopic sample, the photometric redshifts present a rel- 

4 Publicly available at http://code.google.eom/p/EAZY-photoz 
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ative accuracy of Az=(z p h 0 t-Zspec)/(1+Zspec)^3 % with less than 

3 % catastrophic outliers ( i.e. , objects with Az>0.2)._ 

Stellar masses were determined with FAST0 dKriek et al.l 
12009 1 on the U to 4.5 in n PSF-m atched aperture photom etry 
from Pannella et alJ (120141) . using iBruzual & Chariot (120031) de¬ 
layed exponentially declining star-formation histories, which 
mimicks the trend observed with the redshift evolution of 
the SFR-M* main sequence of star-forming galaxies (see 
Pann ella et al.ll2014l for more details). We u se solar metal l icities , 
a Salpeter initial mass function, and the ICalzetti et alJ {2000) 
reddening law with Ay up to 4 magnitudes. 

4. Universal mid-to-far infrared colours of 
star-forming galaxies as a function of redshift 

In order to incrementally take advantage of the better spatial 
resolution of shorter wavelengths from 24 to 500 /im, we start 
by comparing the flux density of individually detected Herschel 
sources to their 24 fi m emission as a function of redshift. The ex¬ 
istence of an IR main sequence regime followed by most galax¬ 
ies, as opposed to a minor population of galaxies in a starburst 
mode, suggests that most distant galaxies do not depart by large 
factors from a typical FIR/24/zm ratio, at any given redshift, 
for each one of the Herschel bands. We therefore extend here 
the IR main sequence concept based on the IR8 = Ljr/L ^ m 
colour index dElbaz et alj|2011[) to the more general relations 
L H erschei/L 24 ^m/(i+z )» where L Hersche i is the luminosity in each one 
of the Herschel bands, i.e. 1. 1 oopm/(i h\ oryvm/ci+z)» ••• 

Our goals are to (i) find the best recipe to construct as realis¬ 
tic mock Herschel images as possible that we will use to test the 
robustness of source extraction against confusion, (ii) determine 
a methodology to benefit from the shorter wavelength Herschel 
bands to improve the source extraction in the longer bands on 
real images, mainly using PACS to deblend SPIRE, (iii) esti¬ 
mate the completeness of this source extraction approach, and 
correct the contribution of faint Herschel sources to the cosmic 
IR background for incompleteness. 

4. 1. Clean galaxies 

The goal of this section is to determine the existence of univer¬ 
sal trends in terms of FIR/MIR colours of galaxies as a func¬ 
tion of redshift. In order to minimise the effects of bright neigh¬ 
bours on galaxies, which would artificially make such relations 
noisier, we have chosen to perform these measurements on a 
sub-population of “clean” galaxies. Clean galaxies are defined 
by exploiting the a priori knowledge of the underlying galaxy 
population on the projected density map, assuming that the lo¬ 
cal spatial distribution of the 24/mi sources also traces the local 
spatial distribution of sources at longer wavelengths. Here we 
define a Clean Index (Cl) as the number of bright 24 /im neigh¬ 
bour sources (S ne i„hbour >0.5S cen , m i) that are closer to the cen¬ 
tral source, S t ± than 12" foj-PACS and 20" for SPIRE. It 
was found by iHwang et al.l 12 010 ) that sources with 2 or more 
“bright” 24yum neighbours systematically exhibit non-physical 
jumps in the FIR SED, resulting from their Herschel fluxes, 
while sources with 0 or 1 neighbour did not. In its original 
definition, the Cl was computed within a distance of 20" by 
reference to the SPIRE 250 /am band. This distance that corre¬ 
sponds to 1.1 xFWHM 25 o Mm being too conservative for PACS, 
we have translated it to the PACS 160/zm band, i.e. 12" (equal 
to 1.1 xFWHM mt , m ). 

5 Publicly available at http://astro.berkeley.edu/~mariska/FAST.html 


4.2. Redshift evolution of the far/mid-infrared colour of 
galaxies 

Here we present a new empirical method based on the correla¬ 
tion of FIR/24 n m colours of individual galaxies as a function 
of redshift to extrapolate the Herschel fluxes to be used to pro¬ 
duce realistic simulations. Libraries of spectral energy distribu¬ 
tions are widely used to fit the FIR emission of galaxies and ex- 
trapolate their bolometric IR luminosity dCharv & Elbazl120011 : 

iDale & Heloull2002t Ida Cunha et al.ll2008h . However, uncertain¬ 

ties in our knowledge about the physics of galaxy evolution and 
not well-constrained measurements in the long wavelengths - 
due to detection and confusion limits- lead to large uncertainties 
when using such libraries, established for local galaxies, to ex¬ 
trapolate FIR fluxes of distant galaxies. After having used SED 
libraries, and found them unsatisfactory to predict the FIR fluxes 
based on 24/j.m as compared to the actual Herschel measure¬ 
ments (in particular the SPIRE data for which extrapolations are 
more extreme in terms of the wavelength gap), we adopt here a 
more phenomenological approach to predict Herschel fluxes by 
fitting the observed evolution of the FIR/24 jim ratio as a func¬ 
tion of redshift, and accounting for the departure of individual 
galaxies from this typical trend based on the closest shortest 
Herschel band. The resulting trends will afterwards be compared 
to standard SED libraries. 

We found that the observed FIR/24 /am colours show well- 
defined trends for each Herschel band as a function of redshift. 
Part of the dispersion of these colours comes from the effect 
of source blending. To minimise this effect in the determina¬ 
tion of the typical colour of individual galaxies, we have split 
our sources in three categories, depending on the presence of 
bright neighbours at 24/mi, within a distance depending on the 
Herschel band considered. This technique has been used to iden¬ 
tify “clean” candidate sources (Hwang et al] |2010t lElbaz et aD 
2010. 1201 ll) as sources with the least probability to be affected 
by the confusion with a bright neighbour, where “bright” desig¬ 
nates sources with a flux density of at least half that of the central 
source. The evolution of the observed SHerschei/S 24pm flux den¬ 
sity ratios as a function of redshift are presented on the left panel 
of Figure Q] In Figure [Q sources defined as “clean” because 
they have zero or one bright 24 gm neighbour are shown with 
red and green dots, respectively. We find that “clean” Herschel 
sources follow a well-defined evolution of the FIR/24 /am colour 
with redshift with an average dispersion of only 0.23+0.07 dex. 
The error bars in Figure Q] represent the median absolute devia¬ 
tion, MAD, i.e. the interval in log io containing 68 % of galaxies 
around the median Herschel/ 24/rm colour. The average of this 
dispersion over the full redshift range is 0.23+0.07 dex, with 
values that vary with wavelength : MAD= [0.17, 0.16, 0.18, 
0.21, 0.31] dex at [100, 160, 250, 350, 500] /im. Note that the 
MAD would have been larger without splitting sources in red¬ 
shift, i.e. MAD around the median colour of all redshifts = [0.20, 
0.18,0.20, 0.33,0.42], 

The sliding median of the sources without bright neighbours 
is shown with a blue line joining blue dots with red centres. In 
order to obtain the purest relations, we used here the filled red 
dots, i.e. the galaxies fulfilling the condition of 0=0 (that is to 
say sources having no bright neighbour), to define the best-fitting 
relations and the dispersion around these relations that are shown 
with the solid blue line in Figure Q] 

The green points with at most one bright neighbour are dis¬ 
tributed symmetrically around this trend and we have tested that 
their sliding median does not present any major departure from 
the one seen with the cleanest sources, i.e. no bright neigh- 


5 





































R. Leiton et al.: Herschel identifies individual galaxies responsible for the 80-290pm CIRB 



10 100 1000 10000 
^24fiin QiJy) 

Fig. 2. Histogram showing the number of 24 gm sources as a 
function of flux density in bins of 0.05 dex (black curve). The 
coloured lines show the number of 24 pm sources in the same 
flux density bins that are detected by Herschel at 100 to 500 gm. 


hour. While the technique allowed us to minimise the impact 
of blending, we are limited to the sources that can be detected 
by Herschel. Therefore, in order to check whether this trend re¬ 
mains reliable to fainter levels, we have compared it to stacked 
measurements. 

4.3. Comparison with stacking measurements 

We have stacked postage stamp images on a selection of 24 pm 
prior positions in order to determine whether the colour trends 
seen for detected sources were also representative of a wider se¬ 
lection or more typical of the brightest sources. For each of the 
sources in our galaxy sample, we produced a cutout, centred on 
the nearest pixel to the sky position of the source, in the PACS 
and SPIRE images of 60x60 pixels, which corresponds to an an¬ 
gular scale of about 10 times the image beam. The cutouts were 
then stacked to create mean images in selected bins of redshift 
and stellar mass. The resulting measurements were corrected 
for flux boosting by clustering using the recipe described in the 
Appendix B2 in lSchreiber et al.l (120141 ). 

We stacked only prior 24 pm positions of sources brighter 
than 50 ply. The reason for this choice comes from the num¬ 
ber of independent beams that can be used within 150 square 
arc minutes with the large SPIRE PSF FWHM. Within this 
area there are about 400, 870, and 1650 beams at 500, 350, 
and 250 pm. respectively. Using all 24pm sources brighter than 
20 piy, meaning about 2400 sources, would imply that we would 
be stacking more positions than there are independent beams in 
the area. Obviously, we would not gain any information by do¬ 
ing this. Instead we decided to stop at 50pjy in order to reduce 
the number of priors to nearly 1000. We note that the Herschel 
detections rapidly drop for 24 pm priors below 200 pjy. Between 
200 and 50 pjy, t hey decrease by a factor 10 for PACS and 100 
for SPIRE at the lElbaz et al.l (1201 ll) detection limits of 1.1 and 
2.7 mJy at 100 and 160pm 5.7, 7.2, and 9 mJy at 250, 350, and 
500pm (see Figure[2}. 

Remarkably, the redshift evolution of the FIR/24 pm colours 
measured using stacking (orange line) follows closely the one of 


the sliding median (blue line and dots), showing that the detected 
sources do not present different SEDs than fainter undetected 
sources probed by stacking down to the limit probed here. 

The narrow dispersion of the clean samples allowed us to 
define relations that best fit the median colours of individual 
sources as a function of redshift for each Herschel band (see ^5] 
for more details). These best-fitting relations are shown with a 
yellow line in Figure Q] The nodes of the yellow lines were cho¬ 
sen arbitrarily to reproduce the sliding median with a minimum 
of relations. 

4.4. Blending and clustering effects 

Sources located at prior positions close to two or more bright 
24 pm neighbours (black open circles in Q} show FIR/24 pm 
ratios more dispersed and systematically higher than the best- 
fitting relations over the whole redshift range. This suggests that 
their colours are strongly affected by blending with nearby bright 
sources, hence are not characteristic of the colours of individual 
galaxies. 

In general, the non-clean sources become confused just by 
the presence of a large projected density of sources at their posi¬ 
tion (by chance associations rather than physical association in 
most cases) and there is no reason for them to have a different 
nature than the clean ones. 

A possible caveat though, would be the possibility that clean 
sources favour sources that are physically less clustered than the 
general population. Indeed, it is known that two effects are play¬ 
ing a role in the confusion of Herschel mosaics: a projected ran¬ 
dom association and a physical association linked to clustering. 
Clustering can boost flux densities by several percents as de¬ 
scribed in the Table B.3 of lSchreiber et alJ d2014t) . The clean in¬ 
dex would preferentially reject those sources and, possibly, bias 
the colours of the resulting galaxies if clustered sources were 
found to be, for example, warmer than non clustered ones. In or¬ 
der to test this effect, we have kept the same number of 24pm 
sources and redistributed them randomly in the field. Physical 
clustering is therefore present only in the real source distribu¬ 
tion and absent in the random one. We do find that the frac¬ 
tion of clean sources rises when distributing randomly the 24pm 
sources. In the real catalogue, 38.8 % of the sources are “clean”, 
i.e. present at most one bright neighbour with >50 % of its 24pm 
flux density within 20", while in comparison this fraction rises 
to 44.6 % in the randomly distributed sample. Concequently, we 
would find 7 % more clean sources in the absence of cluster¬ 
ing. In other words, our clean index has led us to reject 7 % of 
the sources because they were physically associated with bright 
neighbour due to real physical associations. We do not know 
whether such sources are presenting different far-IR colours, but 
their small fraction suggests that their presence in the sample 
would not have strongly modified the colour trends that are used 
here. 


5. An incremental method to predict flux densities 
in the Herschel bands 

We have shown in the previous section that the clean popula¬ 
tion is representative of the whole population, as well as of the 
detected sources since they follow the same trend as the staked 
sources. This also implies that we can use the colours of clean 
detections to constrain a proto-typical IR SED of galaxies or 
at least to extrapolate from 24 pm the flux densities of all prior 
sources in the Herschel bands used here to make the mock im- 
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ages. For comparison, in Figure Q] we also show the main se¬ 
quence SED from lMagdis et al . C2() 1 2) (black dotted line), where 
the correct SED is chosen as a function of redshift, and also the 
main sequence SED from Elbaz et ~afld201 ll) (black dashed line). 
These SEDs globally agree with the colour trends that we find 
here. Noticeable differences appear at z >3, where the SEDs are 
not constrained but sources at such high redshifts make a neg¬ 
ligible contribution to the confusion limit. We note though that 
these large differences are however useful for the identification 
of candidate z >3 sources. 

In the following, we will show that the dispersion of “clean” 
sources around the sliding median of the FIR/MIR relations is 
physical and not due to measurement uncertainties. To demon¬ 
strate this, we will determine the offset of sources with respect to 
the median trend at a given wavelength, e.g. 100pm, then correct 
for this offset the following larger wavelength, i.e. 160 pm. If the 
dispersion is physical and dominated by FIR/MIR differences 
rather than changes in dust temperature, then the dispersion at 
160 pm will be reduced by the same amount. We will make such 
corrections iteratively up to 500pm. After having shown that 
this method does reduce by a substantial amount the dispersion 
of galaxies around the median trends, we will use these median 
trends and the offsets measured at “shorter wavelengths” to re¬ 
versely predict the observed “longer wavelengths” flux densities 
of sources starting from 24 pm to 500pm. This method will al¬ 
low us to attribute realistic observed Herschel flux densities to 
all 24 pm priors and we will then produce mock Herschel images 
from these catalogues that we will compare to the real images to 
verify that they do contain the same amounts of confusion noise 
due to undetected sources. 

5.1. Universal Herschel-24/jm colour evolution and 
dispersion 

The observed dispersion of the colours of galaxies as a function 
of redshift results from the combination of noise in the measure¬ 
ments and physical differences such as PAH equivalent widths, 
distributions of grain sizes, nature (silicates, carbon), and tem¬ 
perature. The upper plot in Figure [3] shows the colour residuals 
with respect to the best-fit relations given in Sect. [3 If noise was 
the dominant cause of this dispersion, we would not expect the 
offset of a source observed at 100 pm to be correlated with the 
one observed at 160 pm. Flux boosting, due to a bright neigh¬ 
bour, would artificially induce such correlation but with a dif¬ 
ferent effect as a function of FWHM and, as seen in Figure Q] 
with a preferential effect of boosting the flux density. However, 
as we will show in the following sections, offsets to the best¬ 
fitting trends are correlated between the Herschel bands sug¬ 
gesting that they are due to physical differences. For example, 
lElbaz et al. ( 201 1) showed that starbursts tend to exhibit stronger 
FIR/MIR ratios than normal star-forming galaxies, which result 
in a boost in all Herschel bands. The effect of dust temperature 
on this boost is a second-order effect. 

Assuming that differences in the FIR/MIR colours of ob¬ 
served galaxies are physical and not dominantly due to noise, 
we can use the e.g. excess of far-IR emission, as probed by the 
100/24 flux density ratio, as an indication that the SPIRE/24 pm 
ratio will also exhibit an excess, hence use this information to 
better predict the larger wavelengths flux densities. 

Since the presence of neighbours results mostly from pro¬ 
jection effects of sources that are not physically connected, the 
selection of clean sources, i.e. sources selected to have at most 
one bright 24pm neighbour, is not expected to bias our analysis. 
Since we allow at most one bright neighbour but reject sources 
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Fig. 3. Residuals of colour index Herschel band minus MIPS- 
24 pm. with respect to the best-fit relations as a function of 
redshift. The upper panel presents the residuals for the ob¬ 
served colour index. The 68 % dispersion around the median is 
equal to [0.17, 0.16, 0.18, 0.21, 0.31] dex for the five Herschel 
bandpasses when accounting for redshift ranges in the best-fit 
equations, while it would have been [0.20, 0.18, 0.20, 0.33, 
0.42] without separating galaxies in redshift slices. Lower panel 
shows the dispersion after accounting for the offset measured 
at shorter wavelengths, e.g. correcting S m^m/S 24 fim by the dis¬ 
tance of S 100 iim/S 24 nm to the median, the dispersion in colour 
indices becomes as small as [0.14,0.12, 0.13, 0.19] at [160,250, 
350, 500] pm, respectively (no correction can be applied to the 
100pm band obviously). The correction process reduces the dis¬ 
persion of the colour distribution because that variation is phys¬ 
ical and not due to noise in the measurements. 


with two or more neighbours, we do not expect our sample to be 
biased against interacting systems of close pairs. 

The resulting best-fitting laws that describe how colours 
change as a function of (1 + z), as shown in Figure Q] are the 
following: 

s wt™ = x 30 -° X (1 + z)- 132 ,z < 1 

= s 24^m X 1 .7 X (1 + z) 2 ' 80 , 1 < z < 1.4 (1) 

= S 2 4t<mXl03x(l+zr L88 ,1-4 < z < 2.3 
= S 24//m X 11.0 ,Z > 2.3 
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where S A and 5 24 ^m are nux densities and S A were 

measured in the Herschel images using 24 gm prior positions as 
in lMagnelli et alJd2009t) using the sample of clean sources. 

These best-fitting laws define the reference colours from 
which any excess or lack of flux of the real Herschel sources 
can be determined: in the former case, a departure from these 
laws can, in principle, mean that the source is gaining flux from 
a companion due to blending (boosting), a sign that sources are 
affected by confusion (alternatively, sources can also be experi¬ 
encing a starburst phase that can resemble a boosted source) or, 
in the latter case, the source may be giving up part of its flux to 
boost others. 

Based on these findings, we believe that we can use the 
S Herschel /S 24 /im best-fit colours to extrapolate Herschel flux den¬ 
sities based on 24 jum measurements for objects below the nom¬ 
inal detection limit of Herschel. 

To construct our list of best-fitted fluxes for all prior sources, 
we calculated the S h f‘ for ~ 2700 sources in the ancil¬ 
lary MIPS catalogue, adopting the spectroscopic redshift where 
available (48%) and otherwise the photometric redshift (37%). 
For the remaining 15% lacking a redshift, we assumed that they 
were faint, distant galaxies with colours following the behaviour 
observed in the high redshift domain in Eqs.|T|to0 


5.2. Iterative correction 

The second step to produce the final flux densities, which will 
be used in the realistic mock images, exploits the advantage of 
the multi-wavelength association of data and the observed colour 
deviations with respect to the best-fitting relations, to iteratively 
correct the colour of clean sources from shorter to longer wave¬ 
lengths. 

We begin by considering the shortest wavelengths that show 
very low or no confusion at all. That is the case for 24 jLim and 
lOOyum sources, at least down to the flux density limit above 
which sources are detected in the 160 /um and SPIRE bands. Here 
we will test the hypothesis that the colour dispersion for clean 
sources around the median colour trend, found in the previous 
section, is physical and not artificially produced by noise and/or 
confusion. If this hypothesis is correct, then we can estimate the 
distance in colour of an individual galaxy to the median colour 
trend at the shortest wavelength, e.g. 100/mi, and use this in¬ 
formation to predict the flux density at the largest wavelengths, 


e.g. 250/rm. In other words, we should be able to reduce the 
dispersion around the median in the SPIRE bands by using the 
colour distance measured in the shorter wavelength bands. 

For any given source detected at 24/jm for which we have a 
redshift determination, we start by determining its S from 
Eq. Q] Then, this extrapolated flux density is compared to the 
observed one and its excess (or deficit) is defined by Cioo as : 


Cioo - S 


be st fit , q. 
lOO/im'^ 


observed 

100/im 


(6) 


where Cioo contains the colour information, from a well sampled 
set of data, that will serve to correct any excess (or lack) of flux 
for that source also detected in a longer wavelength band. In the 
case that the source is not detected at 100/um the constant is set 
to 1, i.e. , we adopt the best fitted 100/rm flux as the flux density 
for that 24 fim prior in the simulations. 

To determine whether offsets at 100 and 160/um are corre¬ 
lated, we define a corrected 160 yum flux density as : 


o corrected _ o observed v 
^ 160/mi “ 13 160/im X Ul0 ° 

_ c observed w ^ best fit t n observed 
~ ° 160/im A ° 100^m /O 100/fm 


(7) 


Figure [TJright presents the distribution of the corrected flux 
densities as a function of redshift. No correction was applied at 
100/im since there is no shorter Herschel band that can be used 
for a correction, hence the upper-left and upper-right panels are 
identical. Below, the 160//m band shows a reduced dispersion 
which confirms that the dispersion seen at 100 and 160/mr are 
not only due to random noise but are correlated. Globally, the 
68 % dispersion around the median of the S ratio is 

equal to 0.16 dex when galaxies are compared to their median in 
redshift intervals, as in Figure Q}left. It becomes as low as 0.12 
dex in Figure [T} right, when correcting for the offsets observed at 
100/rm using Cioo- This implies that one can extrapolate a value 
for S i 60 /jm with a 0.12 dex accuracy with the knowledge of the 
redshift of a source and of its flux densities at 24 and 100/im. 

To extend the process to the following wavelength, 250 fim 
we define Ci6o as the residual offset between the corrected 
160/im flux density and the median best-fitting trend : 


Ci 60 - s 


best fit , q. 
160/im 


corrected 

160/im 


( 8 ) 


Using this new correction factor, one can correct the ob¬ 
served 250/im flux density by combining both Cioo and Ci6o 
as in Eq. [9] and then determine the dispersion of resulting flux 
density distribution around the best-fitting S 250 //m/S 24 /im ratio vs 
redshift trend. The 68 % dispersion around the median gets re¬ 
duced from 0.18 dex to 0.13 dex. We apply a similar procedure 
to correct the 350 and 500 //m flux densities, so that: 


c corrected 
° 250/im 


_ c observed 
~ ^ 250/im 


x Cioo x Cioo 


c observed w r> bestfit t ^ observed 
° 250/im A ° 160/im' J 160/im 


C250 


c best fit j n corrected 
° 250//m ' ° 250/im 


C350 


c best fit / r< corrected 
0 350/im ' 0 350/im 


corrected 

350/im 


= x Cioo x Cioo x C 250 


c observed w cbest fit / r> observed 
° 350/im ° 250/im' ° 250/um 


(9) 


GO) 


( 11 ) 
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c corrected 
° 500/im 


S 500u'm d X C >00 X C 160 X C 250 X C 35 0 


c observed w o bestfit / n observed 
° 500//m ° 350yum' ° 350yum 


( 12 ) 


When no measurement of the source exists in the band k then 
we set the constant as Ck = 1 . 

When applying these iterative extrapolations, we observe 
that the original colour dispersion for clean sources (left panel 
of Figure [TJ is reduced in all Herschel bands (right panel in 
Figure |T}. This means that the colour dispersions are correlated 
from one band to the other, hence the knowledge of the shorter 
wavelengths can be used to improve the extrapolation to the 
largest wavelengths. We can consider those corrected flux den¬ 
sities as good approximations to the real emission of the 24 /jm 
sources in the FIR bands. Figure [3] shows the dispersion of the 
observed and corrected colours for all Herschel bands. 

We find that the use of this iterative offset correction for 
the colours at shorter wavelengths (as in Figure |T} right) reduces 
the MAD of the dispersion of the FIR/24 gm colour with red- 
shift to 0.15+0.04 dex, i.e. [0.14, 0.12, 0.13, 0.19] dex at [160, 
250, 350, 500] fim, respectively (no correction can be applied to 
the 100 yum band obviously), as compared to the initial value of 
0.23+0.07 dex (as in Figure [[[left). Thus, the colour dispersion 
is typically divided by two with our technique, which implies 
that the prediction of the flux density in the most extreme band 
at 500 /urn is known within 55% (within 32% at 250yum), starting 
from the knowledge of the source redshift and flux density in the 
shorter bands. 

To summarise the previous discussion, our aim here consists 
in demonstrating that we can robustly predict the flux density of 
an individual galaxy at wavelength “k +1 ” (e.g. 500yum) on the 
basis of a best-fitting Herschel/ 24/jm colour trend that depends 
on redshift and on a set of colours measured by the coefficients 
Ck (e.g. C ioo, C i 6 o, C 250 . C 350 ) that mimic the variety of IR SEDs 
of individual galaxies. What do not claim to demonstrate that all 
galaxies share the exact same SED, but instead that the measure¬ 
ment of the Ck coefficients can be used to improve the predic¬ 
tion of the flux density of individual galaxies at wavelength k+ 1 , 
where blending is more problematic than in the shorter bands. If 
our assumption was incorrect, we would not observe a reduction 
of the dispersion as we do in Figure [[[right and in Figure [3} 
bottom. The residual dispersion combines the uncertain extrapo¬ 
lation from A < k to A — k with noise due to instrumental effects 
and the source extraction process. In the following section, we 
will use these colours to inject mock Herschel sources in a noise 
map to produce a fully consistent mock Herschel image by ac¬ 
counting for the measured Ck coefficients when building band 
k+ 1 . 


6. Production of realistic mock Herschel images 

6.1. Instrumental noise, PSF, and astrometric uncertainty 

To produce mock images as realistic as possible, all characteris¬ 
tics of the actual GOODS -Herschel maps need to be reproduced 
such as the pixel scale, instrumental noise and point spread func¬ 
tion (PSF) to convolve the input fluxes. 

We start with a map only containing instrumental noise, with 
no source, which we build by dividing the individual scan maps 
of the real PACS and SPIRE observations in two equal sets and 
subtracting one half to the other. To avoid being affected by dif¬ 
ferent shifts in the positions of the half-maps, we combined sep¬ 
arately the odd and even individual scan maps, rather than tak¬ 
ing the first half and the second half which would induce a time 
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Fig. 5. Histogram of the flux distribution of pixels for Herschel 
maps and simulated mock images: observed (black), best-fit 
(red), and iteratively corrected (blue) fluxes. The corresponding 
pixel distribution of the noise maps used in the simulations are 
also presented (dashed lines). The panels below each histogram 
show residuals coming from ratio of the observed Herschel map 
to the iteratively corrected map after subtracting the instrumen¬ 
tal noise as a function of the flux density. The ratio varies around 
1 along the flux bins. The similarity between the observed and 
mock images is remarkable in particular in the flux bins where 
most of galaxies contributing to the cosmic infrared background 
are present (See ^ 8 } 


shift between the two mosaics. The resulting maps are devoid of 
sources, since those were subtracted in the process being equally 
present in both mosaics, but contain the instrumental noise. To 
reach the actual instrumental noise level present in the full mo¬ 
saics, this noise map needs to be divided by a factor twcfl 

Individual sources were then injected by attributing to each 
one of the 2704 MIPS-24/rm priors a flux density in each one 
of the Herschel bands using the method described in Section [5] 
i.e. the Herschel over 24 fim colour at the redshift of the sources. 

6 This factor 2 combines two effects: (i) each half-map contains half 
of the exposure time, hence a noise level that is V2 larger than the real 
noise map; (ii) when we subtract the two half-maps, the noise is com¬ 
bined quadratically hence becomes larger by another factor s/2. The 
two effects combine to make the resulting noise map twice noisier than 
the real map, consequently the signal in that map needs to be divided 
by two to mimic the real instrumental noise map. 
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Fig. 6. Same as Figure [5] but for 250//m (upper panel), 350pm 
(middle pannel), and 500pm (lower panel) 

Sources with no photometric nor spectroscopic redshift were at¬ 
tributed a universal Herschel over 24 pm colour as discussed in 
Section^ 

The mock images were built using the observed PSF of 
the asteroid Vesta for PACS and a Gaussian model for the 
SPIRE PSF, including a random astrometric error of 0.35” to 
the source positions. If we had a perfect knowledge of the po¬ 
sition of the sources and if there was no instrumental noise, the 
confusion noise would be highly reduced by the use of priors in 
the source extraction. This can be tested by generating mock im¬ 


ages with these conditions and sources can then be deconvolved 
with exquisite accuracy. For that reason, even though it is said 
that the confusion comes from the Herschel beam, it is because 
the PSF is not perfectly sampled, its central position linked to 
the prior is uncertain, and that there is instrumental noise on top 
that produce a combined effect that forbids the clear deblend¬ 
ing of nearby sources. In order to include a position uncertainty, 
we have performed a blind source extraction on the real images 
and compared the centroid of the extracted sources to the one 
of the prior positions. This comparison gave an offset of ~0.4” 
rms that we used to randomly offset our mock sources. This un¬ 
certainty results from the combination of the Herschel relative 
pointing error (0.2-0.3”) with the remaining absolute astromet¬ 
ric uncertainty and with the precision of the determination of the 
centres of sources in the Herschel maps. The absolute pointing 
accuracy of Herschel is 0.8” but this was corrected by stacking 
positions of the brightest 24 pm sources on each individual scan 
maps before projecting them in the final mosaic. An offset was 
then computed by comparing the centroid of the detection in the 
stack to the centroid of the postage stamp images, i.e. the centre 
of the actual sources. 

In the following two sections, we demonstrate that the mock 
Herschel images designed from 24 pm priors do contain the 
same amount of confusion noise than the real images and that 
they incorporate more than 95 % of the SPIRE 250 to 500 pm 
sources at the depths reached in the present survey. 

6.2. Test of the realism of the mock images from their pixel 
distributions 

A critical issue in the production of Herschel mock images is 
to make sure that mock sources reach faint enough flux limits 
to reproduce the confusion noise level due to the sea of faint 
undetected sources. Here we have assumed a priori that the list 
of 24 pm priors was large and deep enough to include such faint 
sources making the confusion noise level. In order to test this, 
one needs to compare the histogram of pixel distributions in the 
real and mock Herschel images. 

The close similarity between the real and mock images can 
be seen in Figure [4] Quantitatively, this similarity goes beyond 
significantly detected sources as shown by the comparison of 
the pixel distributions of the real and mock images. The upper 
panels in Figures [5] and [6] show the pixel distribution for the ob¬ 
served Herschel (in black) and mock (red and blue) images. The 
mock images (in blue) were generated using the iteratively cor¬ 
rected flux densities as described in the previous section. Mock 
images using only the best-fitting relations between 24 pm and 
Herschel flux densities as a function of redshift were also built, 
and their resulting histograms are shown (in red). This implies 
that the simulations contain enough faint sources to mimick the 
effect of confusion. Therefore, it is not necessary to include addi¬ 
tional prior sources from other catalogues in shorter wavelengths 
(e.g. optical or UV) to reproduce the background of faint sources 
responsible for the bulk of the confusion noise. The noise intro¬ 
duced in the simulation is also show as a Gaussian distribution 
of pixels (dashed lines). The residual of subtracting the instru¬ 
mental noise to the observed and simulated maps corresponds to 
the light produced by the galaxies in all redshifts in the GOODS- 
North field and its sum is the integrated galaxy light. The lower 
panels in the same figures show the ratio between the noise- 
subtracted distributions of the iteratively corrected image and the 
observed Herschel map. The ratio fluctuates close to a value of 
1 in all the flux bins. The recovery of the this fossil light in the 
simulated and in the real maps is very similar confirming the 
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again the accuracy of the simulation making when the correct 
instrumental noise and sources with realistic flux densities are 
used. 

Moreover, these faint sources are also the ones that make the 
cosmic IR background, hence can be used as a reference to de¬ 
termine the fraction of the background that is actually detected 
as compared to the one present in the images (see ® . This con¬ 
firms the results of Dole et a li ( 2006 ). who found that 24 pm se¬ 
lected sources contribute to more than 70% to the cosmic IR 
background at 10 pm and 160 pm. 

Once the set of realistic mock images are built, it is possible 
to extract their source fluxes and compare them with the input 
data to quantify the level of confusion in each source. Since we 
use the actual position of individual sources, the good quality of 
the mock images indicates that we can use them to not only to 
quantify statistically our source extraction method and the level 
of the confusion noise, but also determine the robustness of flux 
density measurements on individual positions in the images. 

6.3. Are we missing a population of 24pm-dropouts in our 
24pm-based mock Herschel images ? 

A pre-requisite to build the realistic mock Herschel images fol¬ 
lowing the process, as described in Section [6] was to start with 
24 pm seeds. The advantage of this choice is obvious since we 
start from a wavelength probing dust emission with the best spa¬ 
tial resolution, a resolution that is itself based on the IRAC one, 
i.e. at the sub-arcsec level. However, the number of seeds is lim¬ 
ited to 2704 sources detected at 24 pm above the threshold of 
20 ply and, although in theory, it was shown that this threshold 
is enough to detect fainter objects (i.e. less star-forming galax¬ 
ies) than the Herschel bands_at the depths used in the GOODS- 
Herschel survey fsee lElbaz et al. 1201 lh . we now wish to confirm 
that a large population of Herschel, and more particularly SPIRE 
sources, are not absent from our simulated mosaics. 

We already discussed in Section [3] the existence of <1% 
24 run -dropouts that are detected with PACS (tMagdis et al.l 
1201 ll ). We wish here to address the remaining issue of the SPIRE 
sources. Two possibilities are offered to address this issue. As a 
first step, one can study the residual maps after having subtracted 
sources from the SPIRE images at the positions of the 24 pm 
priors. The Herschel residual maps are presented in Figure [7] 
The area used for the catalogue production where the noise level 
is homogeneously and deepest in the PACS bands is identified 
with a blue rectangle in all images. We show for comparison the 
MIPS-24/rm image on the upper-left image. The black spots are 
sources in all images. Note that the SPIRE images cover a much 
larger area than the PACS ones, which explains the presence of 
a large number of sources outside of the blue rectangle. One 
can see that the PACS images contain no clear residual sources, 
while there remains 2-3 spots in the SPIRE residual maps. These 
spots are systematically found close to bright sources in the 
mosaics making it difficult to know for sure whether they are 
truly new sources that would be seen with SPIRE but be 24 /un- 
dropouts, or if they are just residual emission from badly sub¬ 
tracted sources due to the combined effects of noise and uncer¬ 
tain position accuracy of the priors. We therefore decided to gen¬ 
erate another set of simulations to disentangle these two possi¬ 
bilities, i.e. true 24/nn-dropouts vs noise artefacts. 

We designed a more complete simulation that wa s ma de out 
of optical priors following the recipe described in lSchreiber et al.l 
( 2014 ) (see Appendix B) that we briefly summarise hereafter. A 
list of 23672 optical priors were used, an order of magnitude 
more than the 2704 MIPS priors, to generate these simulated 


mosaics. These sources were selected to have a 5-cr detection 
in the near-IR (K$ < 24.5). Stellar masses were derived from 
the fit of their optical to near-IR magnitudes, and then converted 
into a star-formation rate (SFR) using the SFR-M* scaling law 
that has_been observed from z -0 to 3.5 with a 0.3 dex scat¬ 
ter in ISchreiber et al.l (120141) . A small percentage of 2% of the 
sources were attributed a 0.6 dex higher SFR to account for star- 
bursts. We note that we do not use the optical-NIR SEDs to de¬ 
rive SFR since this would imply having a perfect knowledge of 
the extinction and star-formation history of galaxies, and that 
the regions emitting in the optical-NIR are representative of the 
colours of the ones where star formation is taking place. This 
is far from clearly demonstrated and in particular the SFR of 
galaxies, which present an excess with respect to the main se¬ 
quence one, are systematically underestimated. Consequently, in 
order to make sure to have both normal star forming galaxies, as 
well as a population of starbursts that are systematically missing 
from optical-NIR estimates, we decided to distribute SFR ac¬ 
cording to the observed distribution of SFR for each stellar mass 
bin. Then SFRs were converted into a total IR luminosity, Lir(8- 
1000 pm) using the Kennicutf ( 19981) relation SFR=1.72xlO -10 
xL ir . The ICharv & Elba^ t 200llT lR SEDs that were best-fitting 
the stacking data were then used to distribute the Lir to each one 
of the Herschel bands in order to generate a catalogue of input 
fluxes. 

Note that we validated that the sources produced following 
that recipe did also follow the colour trends shown in Figure Q] 
The advantage of this other set of mock Herschel sources is that 
it goes much deeper than the present one. The inconvenient of 
this method is that it only statistically reproduces the Herschel 
sources, i.e. sources are given a SFR that reproduces the dis¬ 
persion of the SFR-M* “main sequence”, while here we try to 
stay as close as possible to the real Herschel fluxes. However 
this deeper simulation can be used to determine the fraction of 
SPIRE sources that are missed above our detection thresholds 
when we limit ourselves to 24 pm priors. This is illustrated in 
Figure [8] where the number of SPIRE sources as a function of 
flux density for each one of the three SPIRE bands are shown 
(black line), as well as the sources that are also detected at 24 pm 
above the threshold of our prior list (red line). As one can see, 
above our detection thresholds of 2.5, 5, and 9 mJy at 250, 350, 
and 500 pm the fraction of sources missed when requiring a 
24 pm prior is at the percent level. 


7. The confusion noise 


Ironically, there is a non negligible confusion around the defini¬ 
tion of confusion noise. It is often defined as the limiting depth 
below which the number of sources gets so large that there is less 
than an arbitrary number of beams per sources. Using this defi¬ 
nition, [Oliver et~ail (l2010h obtain a confusion limit of ~20 mJy 
for the three SPIRE bands defined from the depth correspond¬ 
ing to 40 beams/source in the SPIRE deep surveys and BLAST 
data predicts confusion limits for SPIRE bands of 22, 22, and 18 
mJ y at 250, 350, and 500 pm also equivalent to 40 beams/source 
(iDevlin et al.ll2009l) . 

Alternatively, one can estimate the variance of the sky map 
due to the presence of unresolved sources following the recipe of 
ICondonl d 19741) and la belled as the photom etric confusion noise 


in 


Dole et aiT i 20031) . iNguven et al.l J2010) estimate this noise 


to be as large as 5.8, 6.3, and 6.8 mJy/beam at 250, 350, and 
500 pm, hence leading to lower limits to the 3-cr detection limit 
of 17.4, 18.9, and 20.4 mJy, respectively (lower limits since in¬ 
strumental noise is not included here, but discussed in the pa¬ 
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per). Interesti ngly, th e 40 beams/source definition computed by 
iNguven et al . (12010 1 decreases from 250 to 500 pm going from 
19 mJy at 250 pm to 18mJy at 350 pm, and 15 mJy at 500pm. 


One important caveat in the use of number of beams per 
source is that such systematic definition not only relies on an 
arbitrary number of 40 instead of, e.g ., 20, but more importantly 
assumes that the impact of confusion is the same independently 
of the slope of galaxy counts below this detection threshold. It is 
obvious that if the depth considered is close to the plateau of the 
integral counts, then few sources will appear at fainter flux den¬ 
sities, while if the number of sources below this depth presents 
a steep increase with decreasing flux density, the impact on the 
fluctuations of the sky will be much stronger. 


In practice, what one needs is not a systematic definition 
but rather a quantitative estimate of the photometric accuracy 
that can be reached on most of the extracted sources. This is 
particularly true when confusion can be due to very faint un¬ 
detected sources, as well as relatively bright sources that are 
located at specific positions in the maps, introducing a spatial 
variation of confusion. If we consider the standard Gaussian 
value of an rms of 68 %, we can set our limit to be the depth 
above which more than 68% of the sources (at the faintest 
flux limits) can be extracted with a photometric accuracy better 
than 3 2%, meaning a S/N>3. This definition was introduced by 
ICharv et~ak ( 20041) and later used for both Spit zer and Herschel 
(Ma gnelli et al] 20091 120131 lElbaz et alJl201 fl) . The photomet¬ 
ric accuracy on the extraction of Spitzer MIPS-24pm sources 
was estimated by injecting mock sources on the real maps in 
small number, to avoid overcrowding the field, but a large num¬ 
ber of times. The limiting condition above, which corresponds 
to |A5|/5 ouf = \(S out -S i„)\/S ou t <0.32 for >68 % of the sources, 
was reached at ~20 pjy. At such limiting depth, the number of 
beams per source is as low as 7 beams/source, yet leading to 
reliable photometric measurements for most sources. 


A potential drawback of the latter method is coming from 
the fact that sources are injected at some specific positions 
in the map and even though they can be randomly disposed, 
one cannot estimate the photometric uncertainty on positions in 
the map were bright sources are located. In ultra-deep images, 
bright sources are located all over the image, hence making it 
difficult to study variations of the confusion noise as a func¬ 
tion of position. In order to account for the specific position of 
bright sources in t he Herschel images, iHwang et alJ (120101) and 
Elba z et al.l (1201 Oh noted that there was a correlation between 
the presence of bright sources in the neighbourhood of a specific 
source in the Spitzer 24 pm and Herschel maps which led them 
to select sources using the clean index defined in Section 14.21 
Using the clean index to reject such objects in the GOODS fields, 
lElbaz et akl d201 lb were able to reduce the detection threshold 
accounting for confusion (corresponding to the 3<x limit) by a 
factor ~3 in the SPIRE bands reaching 5.7, 7.2, and 9 mJy at 
250, 350, and 500pm respectively. 


7.1. The Photometric Accuracy Index (PAI) 

For a valid comparison between the confusion levels in the real 
and the simulated maps, not only the mock images have to be 
built as realistically as possible but also the extraction method 
applied must be identical to the method used to produce the 
GOODS -Herschel catalogues. 

The photometric accuracy of the flux extraction in the mock 
images will allow us to quantify with high precision in what de¬ 
gree the extracted sources are affected by bright neighbours and 
by the confusion produced by faint sources. We will use that 
knowledge to explore deeper levels in the detection limits of the 
GOODS -Herschel images. 

To quantify the confusion noise in the simulations, we de¬ 
fine the Photometric Accuracy Index (PAI) as in lMagnelli et al.l 

I2U0I : 

PAI = | AS |/S 

out 1(5 out S in)\/S out • (13) 

In the following, we will consider that sources of a given 
flux density, S can be accurately measured if the extracted flux 
density, S out , of at least 68% of them can be extracted with an 
accuracy better than 32%. We wish to introduce here an impor¬ 
tant distinction between what we will call the “global confusion 
limit’’ and the “local confusion limit”. The global confusion limit 
is the minimum flux density that fulfills the requirement that PAI 
<0.32 for 68% of the sources. We have seen in the previous sec¬ 
tions that we can use the information at shorter wavelengths to 
improve our knowledge at larger wavelengths. In particular, a 
source that is far from any 24 pm bright source has on average a 
Herschel flux density that is less affected by blending. Since our 
mock Herschel sources are located at the same positions as the 
real ones and with similar flux densities, we can use the simu¬ 
lations to reject a fraction of the sources that are more likely to 
be affected by the presence of bright neighbours than others. By 
doing so, we will reduce the effective area of the survey but with 
the advantage of decreasing the flux limit above which output 
fluxes are reliable. In the following, we will call such improved 
flux density limit the “local confusion limit”, local in the sense 
that we account for what happens at the location of sources. 

Figures [9] and ITOlshow AS/S oul as a function of S out for the 
simulations of the GOODS-North field. The plots show how the 
photometric uncertainty increases inversely with the flux den¬ 
sities. Positive values of A S/S oul above 0.32 indicate that the 
extracted fluxes in that flux bin are overestimated. There is a 
characteristic asymmetry in that direction in the SPIRE bands 
that is typical of what one would expect from the contamination 
of bright nearby sources. However, the behaviour of the clean 
sources in these plots will allow us to obtain new depths, pre¬ 
sented below, based on the simulations. 

7.2. Validating the “clean index’’ (Cl) with the “photometric 
accuracy index" (PAI) 


In the following, we will test this method on our realistic 
mock Herschel mosaics, improve it by studying all positions on 
the mosaics and quantify the limiting depth reaching the “3cr 
accuracy for 68% of the sources” at the actual positions of the 
Herschel sources. We will show that one can reach a detection 
limit, that we will from now on call the local confusion limit 
(since it accounts for the local noise due to bright sources as 
determined from shorter wavelengths), 7 times deeper than the 
global confusion limit of ~20 mJy at 250pm (4 times at 350 pm 
and 2.5 times at 500pm). 


In Figures[9]and [10] all sources are shown with a filled black dot, 
while clean sources with at most one bright neighbour (CI<1, 
which is going to be our clean criterion) are represented by open 
green circles. Clean sources, as determined from the 24pm clean 
index (sources with at most one bright neighbour at 24 pm with 
12 and 20"in the PACS and SPIRE bands, respectively), are 
shown with light green dots. Clean sources are nearly indistin¬ 
guishable from the full set of sources in Figure [9] showing the 
100 and 160pm PACS results. On the other hand, the clean sam¬ 
ple tends to be centred on the |A5j/S OMf = 0 line in the SPIRE 
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bands ('Figure UOt . while the rest of the sources, those which hap¬ 
pen to fall close to a bright 24 gm neighbour, are systematically 
above the dotted \ AS\/S 0 ut = 0.32 line. 

This shows that the clean index is a robust method for reject¬ 
ing sources that are heavily affected by confusion. At 250 pm 
for example, one would need to consider a limiting depth of 
~20 mJy to reach the requirement of 68% of the sources with 
\AS\/S a „t < 0.32, which would correspond to the global confu¬ 
sion limit. Instead, we can see that this requirement is reached at 
much fainter levels, down to a depth of 2.5 mJy after rejecting 
non clean sources (vertical dashed line). Using the clean index, 
the local confusion limit is set by the intersection between the 
red envelopes that encompass 68% of the sources around the 
sliding median of the clean sample (blue line) and the horizontal 
dotted lines at |A5 \/S out = ± 0.32, which define our 3-cr limit. In 
these mock Herschel images, we can see that one can expect to 
reach a local confusion limit of 1.1 mJy, 2.2 mJy, 2.5 mJy, 5.0 
mJy, and 9.0 mJy for the 100yum, 160 pm, 250/rm, 350/rm, and 
500 pm Herschel bands, respectively (intersection of the vertical 
lines with the |AS|/S OM , = ± 0.32 in Fig. [9]and ITTTb . These limits 
are slightly deeper than the lElbaz et all tt201 lt> ones at 250 and 
350 pm as discussed in the next section. 


7.3. Comparison of the noise level in the real and mock 
Herschel images 

The PAI can obviously only be determined in the mock Herschel 
images for which we know both S m and S out . However, we can 
check the consistency of the noise level present in the real and 
mock Herschel images by comparing the PAI=\AS\/S ou t t° the 
noise level attributed to individual detections in the real im- 
ages as estimated from the residual maps (see iMagnelli et al.l 
(2013) for a more detailed description). Once detections are ex¬ 
tracted from the real Herschel images, one can estimate the lo¬ 
cal noise at the position of a source by the fluctuations mea¬ 
sured in the residual map. The robustness of the noise estimate 
obtained on residual images has been tested by injecting a lim¬ 
ited number of mock sources on the real Herschel images, and 
extracting their flux density to de rive their |A5|/5 nM> . As dis¬ 
cussed in lElbaz et all (1201 ll) and IMagnelli et al. (1201 3l) . while 
both noise estimators agree for PACS, there is a systematic off¬ 
set for SPIRE. The noise obtained from the residual maps in the 
SPIRE bands is lower than the actual uncertainty in the mea¬ 
surements of the | AS | /S out . Typically, the 5-cr detection thresh¬ 
old estimated from the residual maps corresponds to an effective 
threshold of 3-cr. In order to derive a local noise estimate, one 
needs to use the residual maps and correct the rms by this sys¬ 
tematic offset. 

After accounting for this systematic offset, we can perform 
a comparison of the noise estimated by the PAI in the mock im¬ 
ages (solid red line in Figures [121 131) with the noise estimated 
on the real Herschel images from the residual maps (solid blue 
line). We have here directly applied the strategy described above 
for the real images and, without adding any correction factor, one 
can see that the red and blue solid lines perfectly match, confirm¬ 
ing again that the mock Herschel images are indeed close to the 
real ones. We reach very similar detection limits. In lElbaz et alj 
(20111 the detection limits were based on the blue solid line. 
The red dashed line shows the S/N obtained for the clean sources 
only in the mock Herschel images, we can see that at 250 and 
350/iin the 3-cr limits are slightly deeper (-2.2 and 5 mJy). 


7 . 4 . Representativity of “clean” sources 

We have demonstrated that the clean index (Cl) offered a reli¬ 
able way to pre-select prior positions for which Herschel flux 
density measurements are obtained with a photometric accuracy 
better than 32%. Imposing on the list of 24 pm priors the con¬ 
dition that they have at most one bright neighbour in that band 
is an efficient method to reject galaxies with poor photometry. 
However, one may wonder whether this result is conditioned by 
the fact that we used Herschel images based on 24 pm sources 
only. We showed that the histograms of pixel values of the mock 
and real images were comparable (Figures l9l ITOl). and that the 
S/N ratio as a function of flux density were also comparable 
(Figures 112113b . In this section, we wish to present a comple¬ 
mentary test that reinforces already existing evidence that select¬ 
ing sources by the clean index does not bias their properties, in 
terms of FIR over MIR colour, and to compute the completeness 
of such selection. 

The following test has been carried out for all Herschel bands 
but for brevity is shown here for the 250 pm band alone. We re¬ 
produce in Figure[TT]the 250/24 colour as a function of redshift 
as in Figure Q] The measure of the average 250 pm flux density 
obtained by mean stacking images at the positions of all 24 pm 
sources is the same as in Figure Q] and shown with an orange 
line. It provides a reference for the typical trend of this colour. 
We note that the 250/24 colour presents nearly no variation with 
redshift because, while the redshift increases, the 250 pm band 
goes across the ~ 1 00 pm bump due to the thermalised big dust 
grains, and in parallel the 24 pm band goes across the ~8/rm 
PAH bump. If we select the 24 pm prior positions for which we 
get a PAI-\AS\/S ou t > 0.32 in the mock Herschel image and plot 
in Figure|TT]the 250/24 colour measured on the real Herschel im¬ 
age on those positions, we get a 250/24 ratio, which is typically 
twice as large than the typical one (red solid line in Figure ITO. 
For this test, we have selected the faintest flux density bin, where 
the photometry is more likely to be uncertain, i.e. 3 mJy < S 250 
< 5.7 mJy. Even at such faint flux density levels, we find that if 
we instead select sources that exhibit a PAI-\AS\/S out < 0.32 in 
the mock Herschel image the colour of these sources in the real 
Herschel image is very close to the one obtained from stacking 
as shown by the solid black line in Figure Qj] This shows that 
the simulations do provide an efficient way to identify sources 
with reliable photometry among the sources that are detected in 
the real images. The PAI-\AS\/S oM is a solid reference indicator 
for the local confusion noise but it requires to produce realistic 
Herschel images and one would like to test if a systematic cri¬ 
terion can provide a reliable selection mechanism as well. Clean 
sources, identified with the clean index, are shown with a solid 
blue line with red dots in Figure [IT] They follow with a remark¬ 
able accuracy the trend seen for the PAI-\AS\/S ou r < 0.32 and 
stacked sources. Therefore, selecting sources with the clean in¬ 
dex does not induce a systematic bias in their IR SED. 

The quality of the source extraction done in the mock 
Herschel images is described in the following using two crite¬ 
ria: completeness and purity. The purity defines the photometric 
reliability of the catalogue of sources. It is measured with the 
PAI, i.e. the relative difference between the input flux density 
and the measured one, and sources are considered as “pure” if 
their photometric accuracy is better than 32 %. The dash-dotted 
line in Figure fl4l-tor) shows the purity of the whole sample of 
sources, i.e. using as priors all the 2Apm positions, in the case 
of the 250 pm band. At flux densities lower than 4 mJy, less than 
68 % of the sample is pure. If instead we require that prior posi¬ 
tions be limited to only 24 pm sources with at most 2 neighbours. 
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then this limit is reached at 3 mJy. Sources with CI<1 (i.e. , hav¬ 
ing at most one bright neighbour), are pure for 68 % of them at 
2.5 mJy. Asking 0=0, does not improve this limit on the purity 
of the sample. 

The second criteria, completeness, quantifies the fraction of 
sources detected in comparison with the total number of sources 
that have been injected in the mock Herschel images. The dot¬ 
ted line Figure [T4l-bottom. which presents the CI<2 condition, 
is obviously the one with the largest completeness, but at the 
expense as discussed above of a worse purity. In other words, 
adding sources by allowing to go from 0=1 (dashed line) to 
0=2 (dotted line) enhances the completeness by adding sources 
with a photometric accuracy worse than 32 %, the limit that we 
have here defined as that of “pure galaxies”. The solid red line in 
Figurefl4l-bottom shows the fraction of sources that are detected 
with a PAI<0.32, i.e. , the fraction of pure sources among all 
measurements at the positions of the priors. If the mock Herschel 
images were exactly reproducing the amount of noise that each 
source suffers from then we could use the positions of these 
pure sources to obtain the best completeness in the catalogue 
of sources built out of the real Herschel image. However, we can 
also see that by requiring systematically CI<1 (dashed line), we 
get a completeness that is not so bad in comparison and which is 
much better than the CI=0 request. 

In conclusion, the best compromise to optimise the best pu¬ 
rity for the largest completeness appears to be the requirement 
of CI<1. We note, however, that the completeness correction 
factors estimated from the simulations is only a lower limit, in 
the case of SPIRE particularly, since the mock Herschel images 
were built from a list of 24 gm priors hence miss a potential pop¬ 
ulation of sources below the 24 gm detection limit. 


8. Resolving the cosmic infrared background with 

GOODS-Herschel 

The realistic mock Herschel images have helped us calibrate 
the robustness of the source extraction above a limiting depth, 
where more than 68% of the sources are detected with a pho¬ 
tometric accuracy better than 32%. They also provide a way to 
estimate the completeness of the catalogues extracted above this 
depth as a function of flux density and thus to correct for incom¬ 
pleteness the contribution of Herschel sources to the cosmic IR 
background (CIRB) down to the local confusion limit (~2 mJy), 
which is nearly ten times deeper than the global confusion limit 
(~20 mJy). 


While statistical stu dies provide estimates of the amount of 


EBL due to galaxies (jPatanchon et al 

2009 

Dole et al. 

20061 

Marsden et al. 2009; Bethermin et al. 

2010fc 

Glenn et al. 

2010t 

Bethermin & Dole 2011; Berta et al. 

2011 

; Bethermin et al. 

2012ff they do not allow the study of the properties of the galax- 


ies making this contribution. Using the method of the PA1 and/or 
Cl described above to identify sources at the faintest levels in the 
Herschel maps and correcting their contribution to the EBL from 
their incompleteness as estimated from simulations (as a lower 
limit), we can study a greater population of galaxies making a 
stronger contribution to the EBL than extraction methods relying 
on a single wavelength and at the same time have the possibility 
to study the nature of these individual galaxies. 


8.1. Integrated galaxy light: fraction of the cosmic IR 
background resolved in individual detections 

The combined contribution of individual sources to the CIRB is 
the integrated galaxy light (IGL), which can be expressed in its 
differential form as: 


dIGL dN S band 
dS ~ ~dS X 10 20 


(14) 


where d(V(sr _1 ) is the surface density of sources with flux density 
S b anc i in mJy over a flux bin d.S’ (mJy) and Vband in Hz is the fre¬ 
quency in a given FIR bandpass. The constant 10 20 appears after 
the conversion of 1 mJy = 10 Jl nWm 2 Hz Integrating Ea.[T4l 
down to the Su m , we obtain the IGL in units of nWm 2 sr Our 
IGL measurements correspond to high fractions of the CIRB 
down to very deep flux limits compared with other estimations. 

We have compared the estimates of the amount of extragalac- 
tic background light produced by individual galaxies resolved by 
Herschel coming from the PAI and Cl techniques and found that 
they were consistent within less than 5 %. We present here the 
values obtained using the Cl alone since it does not rely on the 
assumptions on redshifts and SEDs used to produce the simula¬ 
tions. 

We first compute the cumulative contribution of all “clean” 
sources individually detected above the Herschel flux limits - 
numbers in italics on the second line at each wavelength in 
Table[2] Then, we correct these numbers for completeness using 
the correction factors derived from our realistic mock Herschel 
images and obtain the values listed in bold letters in the Table[2] 
For example, over a total of 971 sources injected in the mock 
Herschel image at 250/mi above a flux limit of S J" = 2.5 
mJy, the algorithm detects 398 sources above =2.5 mJy 
that are classified as clean (i.e. have at most one bright neigh¬ 
bour at 24/mi). Therefore, the completeness of our clean list 
is N comp i =41 % in number of sources and IGL comp i =61 % in 
the contribution to the IGL (higher than N comp i since detect the 
brightest sources). 

In comparison with the commonly used definition of the 
confusion limit (e.g. 40 beams/source in iNguven et al.l l201( )i). 
we have also computed the equivalent of our confusion lim¬ 
its in numbers of beams/source. Various papers have used dif¬ 
ferent definitions of the beam size, here we use a_definition of 
the beam ar e a mat ching the ones used in iDole et all (l2003h and 
lOliver et~ai] ((2010) for consistency in the comparison: Abeam = 
7 r x (0.73xFWHM) 2 . The beam size at 250/mi is 0.154 arcmin 2 
and the number density of clean detections is equivalent to 3 
beams/source, knowing that this level was reachable only by the 
use of the shorter bands to allow the identification of robustly 
determined photometric measurements. We list in Table [3] these 
numbers for all 5 Herschel bands. 

The resulting background estimated in nW m _2 sr 1 are then 
compared to estimates of the extragalactic background light 
(EBL) at these wavelengths, in order to determine the resolved 
fraction of the CIRB. Two options exist for this comparison 
that consists in using as reference values either the COBE mea¬ 
surements or statistical extrapolations of Herschel counts using 
stacking and P(D) analysis. Error bars on the COBE values are 
larger, especially at 100/mi, where the correction for zodiacal 
light is hi ghly uncer t ain, he nce we use here extrapolated values 
from ( a ) iBerta et all (1201 ll) and by ( b ) iBethermin et alJ (12012 1) 
(listed in the last column of Table [2]) as a reference to compute 
fraction of the resolved CIRB. 

InthePACS 100 and 160//m bands, sources brighter than 1.1 
and 2.2 mJy make around 50-60% of the CIRB at these wave- 
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wavelength 

(yum) 

S > S 
(mJy) 

resolved IGL 
(nWm~ 2 sr 1 ) 

fraction of CIRB 

Comments 

reference CIRB* 
(nWm 2 sr') 

100 

1.1 

6.9 

54% 

after completeness corr. 

12.61 [-1.74;+8.31] n 

100 

1.1 

6.1 

48% 

before completeness corr. 


100 

1.2 

8.35 

66% 

Berta et al. (2011) 


100 

0.54 

9.45 

75% 

Magnelli et al. (2013) 


160 

2.2 

8.2 

60% 

after completeness corr. 

13.63 [-0.85;+3.53]" 

160 

2.2 

7.2 

53% 

before completeness corr. 


160 

2.0 

9.49 

70% 

Berta et al. (2011) 


160 

1.29 

10.27 

75% 

Magnelli et al. (2013) 


250 

2.5 

6.0 

55% 

after completeness corr. 

10.13 [-2.33;+2.60] 6 

250 

2.5 

3.4 

34% 

before completeness corr. 


250 

20 

1.73 

17% 

Oliver et al. (2010) 

[resolved sources] 

250 

20 

1.55 

15% 

Bethermin et al. (2012) 

[resolved sources] 

250 

2.0 

7.4 

73% 

Bethermin et al. (2012) 

[using stacking] 

350 

5.0 

2.1 

33% 

after completeness corr. 

6.46 [-1.57;+1.74]* 

350 

5.0 

1.5 

24% 

before completeness corr. 


350 

20 

0.63 

10% 

Oliver et al. (2010) 

[resolved sources] 

350 

20 

0.77 

12% 

Bethermin et al. (2012) 

[resolved sources] 

350 

2.0 

4.50 

70% 

Bethermin et al. (2012) 

[using stacking] 

500 

9.0 

0.37 

13% 

after completeness corr. 

2.80 [-0.81;+0.93] 6 

500 

9.0 

0.30 

11% 

before completeness corr. 


500 

20 

0.15 

5% 

Oliver et al. (2010) 

[resolved sources] 

500 

20 

0.14 

5% 

Bethermin et al. (2012) 

[resolved sources] 

500 

2.0 

1.54 

55% 

Bethermin et al. (2012) 

[using stacking] 


Table 2. Summary of the resolved IGL, its contribution to the CIRB down to specific flux density limit S u m in each SPIRE band for 
the GOODS-North field. The percentages correspond to the fractions of the IGL with respect to the total extrapolated estimations 
of the CI RB reported by ( a ) lBerta et all (1201 ll) and by ( b ) iBethermin et all (1201 2l) . The new flux limits S u m in the PACS bands are 
similar to Elba z et ail (1201 ill but the SPIRE depths are deeper. 


wavelength 

(pm) 

(1) 

5 lim 
(mJy) 
(2) 

Abeam 

(arcmin 2 ) 
(3) 

Ncompi 

(%) 

(4) 

IGL comp i 

(%) 

(5) 

Nreai 

(6) 

beams/source 

(Nbeams!Nreal) 

(7) 

100 

1.1 

0.021 

80 

89 

568 

12 

160 

2.2 

0.057 

78 

88 

589 

4 

250 

2.5 

0.154 

41 

61 

342 

3 

350 

5.0 

0.292 

58 

72 

185 

3 

500 

9.0 

0.632 

78 

81 

43 

5 


Table 3. Summary of the number of sources and beams involved in the simulations, as well as in the Herschel maps. Col.(2): Flux 
density above which sources are detected in the simulation, i.e. confusion limit defined as the depth above which 68 % of the sources 
get measured with a photometric accuracy better than 3 2%, Col.(3): Abeam is the beam area at the central wavelength listed in Col.C 1) 
defined as n x (0.73xFWHM) 2 (iDole et al,ll2003fc [Oliver et al.ll2010l) . Col.(4): N comp i is the completeness in % in number of clean 
sources (neib2A <1), i.e. the % of sources injected with S m >Su m that are detected in the simulation as clean and with S 0Ut >Su m . 
Col.(5): IGL comp i is the completeness in % in contribution to the IGL. Col.(6): N rea i is the number of clean sources detected in the 
real Herschel maps above S n m . Col.(7): beams/source is the number of beams per source, where Nbeams= 150' 2 /Abeam- 


wavelength 

(/mi) 

S > S n m 
(mJy) 

(z) 

( Z)f 

<M*> <M*)/ 

(10 u M o ) 

100 

i.i 

0.94 

0.86 

0.8 

1.4 

160 

2.2 

0.96 

0.98 

0.8 

1.6 

250 

2.5 

0.96 

1.10 

0.9 

2.0 

350 

5.0 

0.94 

1.17 

1.2 

2.4 

500 

9.0 

1.26 

1.54 

2.3 

3.6 


Table 4. Summary of the median (z) and flux-weighted average redshift (z)f, as well as the median (M*) and the flux-weighted 
stellar mass (M*)/ of the galaxy population that produce the resolved IGL down to specific flux density limit S /„„ in each SPIRE 
band found based on the simulations of the GOODS-North field. 


lengths as compared to the ~75 % obtained by iMagnelli et al.l 
( 2013 1 counting sources down to 0.6 and 1.3 mJy in GOODS- 
South. However, thanks to our Herschel mock image we are able 
to correct these values obtained for clean sources for incomplete¬ 
ness. After having applied this correction, we get fractions of 


the resolved CIRB of 54 % and 60 %, which are comparable to 
the ones obtained using stacking and P(D) analysis as shown 
in Figure [15] This is also true for the SPIRE bands, since in 
both cases, PACS and SPIRE, our estimates shown with solid 
coloured lines (after completeness correction) can easily be ex- 


15 



































R. Leiton et al.: Herschel identifies individual galaxies responsible for the 80-290pm CIRB 


tended by eye to the symbols showing the locus of statistical 
estimates using stacking and P(D ) analysis. 

The most remarkable improvement of our method can be 
seen in the SPIRE bands, where we resolve 55, 33, and 13% 
of the background light at 250, 350, and 500 fim with individ¬ 
ually detected sources thanks to the use of our estimate of the 
“local con fusion noise” as co mpar ed to 15, 12, and 5 % as mea¬ 
sured by lOliver et~akl (1201 Ol) and IBethermin et al.1 (1201 2l) with 
the global confusion limit of ~20 mJy. Our approach has allowed 
us to resolve in individual galaxies 3.5 times more extragalactic 
background light at 250 fim and more than ~3 times at 350 and 
500 yum. 

We present in Fig. [15] our estimates of the IGL both before 
(dashed lines) and after (solid li nes) correcting for complete¬ 
ness. D i rect measurements fr om iBerta et ail d2010^ Berta et aTI 
d20nh. dMagnelli et al.ll 20 Tl forPACS and (lOliver et al.ll2010l) . 
Bethermi n et al.l(l201 2l) for SPIRE are show for comparison with 
filled symbols, while extrapolations with stacking and P( D) 
analysis are shown with open symbols. We also present expec¬ 
tations of the IGL fractions predicted for Herschel based on 
a model of multi-wavelenth galaxy counts by Le Borane et al. 
120091) . 


We compared our resolved IGL fractions to the ones pre¬ 
dicted by the backward evolution model of iBethermin et al.l 
( 201 i ll), which reproduces galaxy counts in all Herschel bands. 
In this paper, the authors calculate the confusion limit expected 
for different telescope diameters in the FIR using the criterion 
defined in lDole et al.l (120031) by requiring sources to be brighter 
than 5-crP hotometric , where o-p 11010 ™ 1 ™ j s the noise due to the fluc¬ 
tuations induced by the sources below this detection limit, and 
such that less than 10 % of the sources have a neighbour closer 
than 0.8xFWHM above that detection limit. The standard global 
confusion limit fits the predictions of the model for a 3.5m class 
telescope, but our improved depths in the 250 pm band corre¬ 
spond to the standard confusion limit of a telescope with an aper¬ 
ture close to ~10m (see Figure 14 and Tab. 8 in iBethermin et al.l 

Hon]). 


8.2. Nature of the sources making the CIRB 

We present here the characteristic redshift and stellar mass of 
the galaxies that make the resolved fraction of the CIRB de¬ 
scribed in the previous section. We provide two methods to es¬ 
timate these characteristic values, the median - (z) and (M*) - 
and flux weighted average - (z) f and (;V/*)y computed following 
Eas. fl5ll6l 


<z)f = 


E/-S/U) xzi 
X,-S\(T) 


(15) 


peak of the IR background around 140-160/rm is expectedly 
produced by galaxies around z~0.9. Increasing the wavelength, 
sources at larger redshifts contribute a larger contribution with 
respect to this redshift, hence we observe the increasing redshift 
to 7~ 1.3 in the higher wavelength bands. 

Our median redshifts in 100 and 160/tm are consistent with 
the results of Magnelli et al.l (1201 3l) using the PEP/GOODS- 
Herschel data who obtained z = 0 . 85^033 for sources S loo/mi 


>1.5 mJy, and z=0.94 + np for sources S m/mi > 2.5 mJy (see 
also [Berta et al. 2010, 2011) . In the case o f the SPIRE bands 
(see Figure [161-left). IBethermin et al.l (120121) found that half of 
the CIRB is emitted at the typical redshifts of 1.04, 1.20, and 
1.25 for 250, 350, and 500 yum, respectively, by staking data in 
COSMOS and GOODS-north as part of the HerMES. This in¬ 
creasing trend in redshift in the higher wavelengths band is also 
consistent with our results in the SPIRE bands, while we find 
larger median redshift in 500 yum. 

Thanks to our method, it is now possible to study the in¬ 
dividual galaxies responsible for the peak of the CIRB around 
250 yum. The major difference and strength of our approach with 
respect to other methods (e.g. stacking) is the fact that we are 
identifying the contribution of each individual galaxy to the 
background instead of studying the average contribution of pop¬ 
ulations of galaxies. This is im portant for studies bas ed on large 
galaxy samples. For example, Daddi et aD (i2007bllal) reported 
that, based on stacking, z =2 sources with a “mid-IR excess” 
were powered by AGNs in the mid-IR. However, this excess was 
coming only from a population of about 30% of the sources that 
were indeed powered by AGNs and numerous enough to make a 
strong signal in the stacking. It was impossible from the stack¬ 
ing alone to know if the characteristics found could be applied 
to 30% or 100% of the sources. 

We present in Figure [16] their redshift and stellar mass dis¬ 
tributions. Both distributions are well peaked, and it is quite in¬ 
teresting to see that the stellar mass distribution peaks at a value 
of M* ~9x 10 10 M 0 , close to the mass of the Milky Way with 
a narrow dispersion of only 0.3 dex in logio- Interestingly the 
same happens in the other bands, which present a similar char¬ 
acteristic mass. The redshift peak is at z ~0.95. 

In consequence, most of the cosmic infrared background is 
due to the radiation of siblings of the Milky Way around 7.5 
billion years ago. This finding is in agreement with the fact 
that the stellar mass function of star-forming galaxies does not 
seem to vary strongly with redshift up to z ~ 2 dBell et al.ll2003l 
12007b iPozzetti et al.ll2010[ lllbert et akl 1201 Oh and with the sta¬ 
tistical anal ysis done using stacking of galaxy populations by 


tistical anal ysis dor 
lYiero et akl (120 13l) . 


9. Conclusions 


<M*)y = 


ZiSiWxM+j 

ZiSiU) 


(16) 


where S are the flux densities of the sources contributing to 
the CIRB for a Herschel band centred at A. 

The IGL in the PACS and SPIRE bands to the CIRB, as well 
as the typical redshifts and stellar mass of the galaxies contribut¬ 
ing to the CIRB are summarised in Table [2] We see that the typ¬ 
ical redshift of the galaxies responsible for a dominant contri¬ 
bution to the infrared background light increases with increasing 
wavelength. This is naturally expected as a result of ^-correction, 
since the peak emission of the IR spectral energy distribution of 
individual galaxies, like e.g. the closest starburst galaxy M82, 
lies around 80 fim (with a moderate 7)/ uv , variation). Thus, the 


The main conclusions of this work are summarised below: 

- We have developed a method to robustly predict FIR fluxes 
based on the knowledge of their redshift and emission at shorter 
wavelengths starting from 24 /urn. 

- Using those predicted FIR fluxes, we were able to build 
realistic simulated Herschel images of the GOODS-North that 
we showed to robustly reproduce the confusion noise present in 
the real Herschel maps. 

- We showed that the extrapolation of the FIR emission of 
24 fim sources detected down to 20 fiiy make up enough signal 
in the mock images to sharply reproduce the pixel distribution of 
the observed Herschel maps, suggesting that the bulk of the con¬ 
fusion noise, even fainter than the detection limits, is included in 
the mock images. 
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- We used the difference between the extracted photometry 

of sources in the simulated images and the input fluxes used to 
make the mock images - the Photometric Accuracy Index ( PAI ) 
- to determine the degree of local confusion in the GOODS- 
Herschel maps, i.e. the photometric uncertainty at each specific 
position on the mosaics. _ 

_ - We validate the C lean Index (Cl) approach (Hwang et alJ 

l20HlElbaz et ali201(H) that identifies clean and confusion dom¬ 
inated Herschel sources from the local environment of 24 jim 
sources clean sources in the GOODS -Herschel data based on the 
quantification of the number of contaminating close neighbours 
to a given source. We confirm that the optimal selection rule for 
clean sources is when CI< 1. While the Cl identifies sources that 
respect the PAI criterion (|AS|/S oa/ <0.32), the mock Herschel 
images combined with the PAI allow one to identify reliable 
sources with a stronger completeness and to set deeper flux den¬ 
sity limits for the data in the SPIRE bands. 

- After correcting for incompleteness, the emission of 
sources detected with our approach in the GOODS-North field, 
we find that galaxies individually detected by Herschel down to 
1.1 mJy at 100;um and 2.2 mJy at 160//m make up as much 
as 54% and 60% of the extragalactic background light, whereas 
down to 2.5, 5, and 9 mJy sources produce 55, 33, and 13.0% of 
the EBL at 250 //m, 350 fi m, and 500 gm, respectively. 

- Interestingly we find that the dominant contributors to the 
CIRB in all Herschel bands appear to be galaxies analogous 
to the Milky Way in stellar mass with M* ~9x 10 10 M G rein¬ 
forcing the idea of a dominant dark matter halo mass in which 
stars formed over the Hubble time (IConrov & Wechsiei1 f2009: 
iBethermin et al.ll2013h . 
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Fig. 9. Photometric Accuracy index ( PAI ) , A S /S ou , = (S out - 
Sj„)/S out, as a function of the output flux S out . The median trend 
of |AS \/S out is shown by the blue line and the red envelope lines 
contains the 68% of the clean galaxies, i.e. sources associated - 
with a 24pm source having a CI< 1 (green dots). Dotted horizon- 
tal lines at |AS \/S out = ± 0.32 mark the limits that define sources !§ 
having good photometry. The intersection between the 68% en¬ 
velope line and the good photometry limit defines the depths for 
Herschel data based on the simulations (vertical dashed lines; 

1.1 mJy, 2.2 mJy, 2.5 mJy, 5.0 mJy, and 9.0 mJy for the 100pm, 
160pm, 250pm, 350pm, and 500pm Herschel bands, respec¬ 
tively). The histograms on the bottom and right sides of each 
panel show the distribution of the whole sample of galaxies 
(dotted line), the clean galaxies (solid green line) and the com¬ 
plementary non-clean galaxies (solid black line) as a function 
of the flux density (horizontal histogram), and as a function 
of AS/S out (vertical histogram). The comparison of the solid 
green and dotted lines shows that the proportion of clean sources 
remain identical in all flux density bins, i.e. , there is no bias 
against faint or bright sources when dealing with clean sources 
alone. 
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Fig. 10. Same as Figure |9]but for 250pm (upper panel), 350pm 
(middle pannel), and 500 pm (lower panel) 
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Fig. 11. Comparison of observed S250/S 24 colours for sources 
with flux densities between 3 < S250 < 5.7 mJy selected by 
their associated photometry accuracy from the simulations. The 
median colours of sources with a good photometric accuracy 
|AS|/5 oa , < 0.32 (black line) follow the colour trend, of clean- 
selected sources Cl < 1 (blue line) and the stacked sources (or¬ 
ange line) while the median colours of sources having poor pho¬ 
tometric accuracy |A5|/5 oa , > 0.32 (red line) are systematically 
displaced to brighter 250//m fluxes, an indication of blending. 




Fig. 13. Same as Fig.[12]but for the SPIRE bands. 
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Fig. 12. Signal-to-noise ratio as a function of flux for extracted 
sources in the HerschelPACS maps (blue), the simulated images 
(solid red) and for the subset of clean sources in the simulation 
(dashed red line). The horizontal lines marks the 3-cr level. The 
similarity in the trends demonstrates that the simulations repro¬ 
duce well the sources of noise in the real image, i.e. the faint 
population of 24 /im prior sources, as well as the instrumental 
noise. Flux limits can also be defined for the Herschel data by 
the 3-cr limits reached by clean sources (intersection between the 
dashed red line and the horizontal line). 
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Fig. 14. Plots presenting the quality of the source extraction on 
the mock Herschel images from a “purity” (top) and “complete¬ 
ness” (bottom) criterion. We present here the example of the 
250 pm band, top: purity of the sources extracted in the mock 
images defined as the fraction of sources with accurate photom¬ 
etry, i.e. A.S' \/S, mt < 0.32, among the sources selected using the 
clean index, CL The different lines present the sources with no 
bright 24 pm neighbour, CI=0 (solid line), at most one 24 pm 
neighbour, CI=1 (dashed line), at most 2 neighbours, 0=2 (dot¬ 
ted line), and the full sample of sources without any filtering 
(dash-dotted line), bottom: completeness of the samples selected 
using the clean index (same definition of the lines), defined as 
the ratio in percentage of the number of sources filtered by the 
clean index over the total number of sources in each flux den¬ 
sity bin (differential completeness). The red solid line shows the 
completeness of the sources with PAI=\AS\/S out < 0.32, which 
are by definition associated with a 100% purity. 
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Fig. 1. Redshift evolution of the mid-to-far infrared colours of Herschel galaxies. Each plot shows the ratio of the observed Herschel 
over MlPS-24//m flux densities as a function of redshift. Figures on the left column show the colours from the observed Herschel 
measurements. Figures on the right column show the same dataset but after having corrected the 160 to 500 /im over 24/rm ratios 
of the individual galaxies by the offset observed in the 100 to 350 i-im shorter wavelength bands. The points indicate the positions of 
individually detected sources separated in 3 populations, depending on the presence of bright 24/j.m neighbours (i.e. with a S 24 /m , 
greater than half of the central source) in their close environment (at a distance of 12"for PACS and 20"for SPIRE). This number is 
the “clean index’’. Cl. Filled red dots: galaxies with no bright neighbours (0=0); filled green dots: galaxies with at most one bright 
neighbour (CI<1); open circles: remaining 3-cr sources (CI>2). The sliding median of the 0=0 galaxies is shown with a blue line 
and blue circles with central red dots. The trends for the mean stacking at the positions of the 24/jm prior sources brighter than 
50/iJy are illustrated with orange lines. The yellow line corresponds to the best-fit relations described in kq. |T| to [5] of §|5] that fit 
the sliding median. The colours obtained for the I R spe ctral energy distributions of the library of lMagdis et al ( 2012 1 (black dotted 
line) and the main sequence SED of lElbaz et al.l ( 2011 1 are shown for comparison.) 
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Fig. 4. Selected area of the real (left) and the simulated maps (right) for GOODS-North field in 250pm. Each stamp is ~9’x9’ (see 
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Fig. 7. Residual images in GOODS-North field after source extraction using 24pm priors in each band. The maps show how efficient 
is the extraction process in the area covered by the prior catalogue (blue rectangle), showing only very few cases where the remnants 
of bright sources badly subtracted are still present in the SPIRE residual maps (black spots). 
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Fig. 8. Number (lower panels) and fraction (upper panels) of 24 gm sources having a SPIRE counterpart (red) and those with no 
associated SPIRE detection as a function of the flux density in 250 yum (left panel), 350/Jin (central pannel), and 500 /urn bands (right 
panel). Vertical lines show the flux limits reached using the simulations found in the present work 2.5 mJy for 250/rm 5.0 mJy for 
350/rm, and 9.0 mJy for 500gm. (see § 17. II for more details of how these flux limits were obtained). Consequently, the missing 
24 /urn population at the flux limits of the simulations is very marginal. 
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Fig. 15. Cumulative contribution to the CIRB of the IGL as function of the flux density limit in 1 00 /mi (violet), 1 60 /mi (blue), 
250/rm (green), 350pm (orange), and 500pm (red) in the GOODS north field. Dashed lines show the cumulative flux density 
directly from the observed GOO DS -Herschel sources and the solid lines are the same but corrected by co mpletenes s from the 

erta et al.l201 lllMa gnelli et alJ2013h . stacked (iBerta et al.1 


simulations in this work. Res olved (IBerta et al. 201(1: Oliver et al.l201 (i lBi 

1201 Ol iBithermin et alll2012h and P(D) results (IBerta et all 201 11) at different adopted S n m are also shown for comparison, as well 
as IGL predictions for Herschel based on a model of multi-wavelength gal axy counts (ILe Borgne et alj|20 09). Horizontal dashed 
bands show the total extrapolated CIRB in the SPIRE bands adopted from lBethermin et all ( 2012 ) using the same colour coding 
than for the curves and their widths represent the limits of the errors. In the case of PACS total values, vertical thick bars mark the 
error in the extrapolations of the CIRB calculated bv lBerta et al.l d20 1 ll) . 
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Fig. 16. Distributions in redshift (left) and stellar mass (right) of the sources responsible for 55 % of the CIRB at 250 yum. Vertical 
dashed lines show the median of the distributions. 
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